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ABSTRACT 


A 


The  study  describes  a number  of  algorithms  for 
solving  single  facility  deterministic  location  problems 
in  which  the  planar  assumption  is  not  appropriate.  Trans- 
formations on  the  non-Euclidean  spherical  space  are  com- 
bined with  efficient  solution  techniques  in  En . Extensive 
use  is  made  of  projective,  synthetic  and  analytic  geometry. 

Two  algorithms  are  presented  for  solving  single 
facility  problems  with  the  objective  of  minimizing  the 
total  sum  of  costs  (minisum) . Due  to  the  non-convex 
nature  of  the  problem,  a local  optimum  is  obtained.  Com- 
putational experience  in  solving  a number  of  test  problems 
is  reported.  Theoretical  results  concerning  narrowing  of 
the  search  region  are  presented  as  well  as  a number  of 
special  properties  of  the  problem.  Application  of  the 
single  facility  results  to  the  location-allocation  class 
of  multifacility  minisum  problem  is  investigated. 

Algorithms  for  the  single  facility  problem  having 
the  objective  of  minimizing  the  maximum  distance  (minimax) 
are  also  considered.  It  is  seen  that  existing  algorithms 
in  E^  can  be  used  to  solve  the  problem  when  all  demand 
points  are  containable  in  a hemisphere.  Application  of 
single  facility  results  to  solution  of  a specific  type  of 
multifacility  minimax  problem  is  investigated. 
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CHAPTER  I 
INTRODUCTION 

I . 1 Background  of  the  Research 

Human  concern  about  distance  is  an  integral  part  of 
everyday  life.  Frequently  decisions  for  or  against  an 
alternative  are  primarily  based  upon  the  constraining  ele- 
ments of  "too  near"  or  "too  far."  The  theory  of  optimal 
facilities  location  provides  tools  which  can  assist  one 
in  selecting  sites  for  servicing  facilities  that  will 
directly  affect  the  distance  traveled. 

Cooper  (1963)  credits  Cavalieri  with  first  con- 
sidering, in  1647,  the  problem  of  finding  a point  at  which  the 
sum  of  those  distances  from  three  given  points  is  a minimum. 
This  is  referred  to  variously  as  the  Steiner  or  Fermat 
problem  in  the  literature.  Sylvester  (1857)  evidently 
first  proposed  the  problem  of  finding  the  points  in  the 
plane  such  that  the  maximum  distance  to  a finite  set  of 
points  in  the  plane  is  minimized.  Although  location  theory 
has  developed  far  beyond  these  fundamental  mathematical 
problems,  new  observations  and  results  are  still  being 
reported  concerning  them  (Sokolowsky  1976) . 
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A new  era  in  location  theory  began  with  development 
of  an  iterative  process  to  solve  a generalization  of  the 
Steiner-Fermat  problem.  The  technique,  developed  inde- 
pendently by  Cooper  (1963)  and  Kuhn  and  Kuenne  (1962), 
was  first  proposed  by  Weiszfeld  (1936) . In  the  past  twenty 
years  literally  hundreds  of  papers  have  appeared  on  the 
subject  in  many  different  fields.  Regarding  the  general 
literature  of  location  theory  as  it  exists  today,  ReVelle 
et  al.  (1970)  found  it  convenient  to  consider  two  main 
structural  categories: 

1.  Location  on  a plane,  characterized  by  an  infi- 
nite solution  space  and  a distance  measurement  according 
to  any  particular  metric. 

2.  Location  on  a network,  characterized  by  a solu- 
tion space  consisting  of  points  on  a network  and  the  dis- 
tance or  time  measurement  d^ . is  the  length/time  of  the 
shortest  path  from  node  i to  node  j . 

This  morphology  is  generally  accepted  as  repre- 
sentative of  location  systems  theory. 

Early  maps  of  the  world  (circa  eleventh  century) 
depicted  the  earth  as  a flat  plane.  The  discovery  that  the 
earth  was  essentially  spherical  created  problems  for 
cartographers  that  are  still  being  investigated  today. 

It  is  well  known  (Hilbert  1952)  that  there  exists  no  iso- 
metric (length-preserving)  transformation  from  a sphere 
to  the  plane.  This  property,  or  lack  thereof,  eliminates 
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the  possibility  of  a transformation  so  that  Euclidean 
geometry  can  be  used,  in  the  large,  to  measure  spherical 
dis  tances . 

Due  to  the  relative  intractability  of  working  in 
non-Euclidean  space,  however,  much  effort  has  been  expended 
in  the  past  in  attempting  to  develop  transformations  from 
non-Euclidean  to  Euclidean  space  which  preserve  desired 
properties.  A review  of  such  attempts  is  provided  by  Angel 
and  Hyman  (1972).  Also,  in  another  vein,  research  is 
ongoing  in  attempts  to  create  optimal  transformations  from 
the  sphere  to  the  plane  which  minimize  the  error  of  a given 
parameter  (area,  distance,  etc.)  for  a specific  region. 
Recent  work  has  been  done  by  Milnor  (1969)  and  Gilbert 
(1974) . 

In  spite  of  the  foregoing  efforts  there  remain  prob- 
lems in  location  theory  for  which  the  Euclidean  assumption 
is  totally  inappropriate.  Unfortunately,  though,  the 
majority  of  theoretical  developments  in  several  related 
fields  (location  theory,  theoretical  geography,  and  regional 
science)  require  a geometrical  framework  based  upon  the 
assumption  of  a Euclidean  plane.  Considering  the  nature 
of  the  earth,  it  is  surprising  that  so  little  considera- 
tion has  been  given  to  non-Euclidean  bases,  specifically 
the  great  circle  metric.  It  is  such  a situation  to  which 
this  research  is  addressed. 

A simple  illustration  of  the  magnitude  of  error 
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which  can  occur  is  provided  in  Table  1.1  for  a three  point 
Steiner  problem  on  the  sphere  (see  Data  Set  D9  of  Appendix 
B) . Using  both  a Euclidean  assumption  and  a spherical 
assumption,  the  problem  was  solved  by  an  existing  algorithm 
for  the  plane  and  methods  developed  in  this  thesis  for  the 
sphere.  The  coordinates  of  the  points  are  in  degrees 
latitude  and  longitude.  For  the  planar  problem,  the 
spherical  coordinates  were  used  as  Cartesian  coordinates. 


I 
I 

When  using  the  planar  assumption  there  is  an 
18.85°  error  in  location  of  the  optimum  facility,  or  approx- 
imately 1300  miles  on  the  earth's  surface.  Also,  there  is 
over  a 307o  error  in  the  optimum  objective  function  value. 

Due  to  the  increasingly  global  nature  of  corpora- 
tions and  international  defense  structures,  it  would  seem 
useful  to  consider  large  region  location  problems.  As 
Warntz  (1966)  declares,  perhaps  the  time  of  "community 
earth"  is  not  too  far  removed  into  the  future.  Due  to  the 
infeasibility  of  developing  a transformation  to  handle 

! - 


Location 

Obiective  Function  Value 

(Lat . , Long . ) 

Great  Circle  Metric 

Euclidean  Norm 

Optimum 
on  Sphere 

(80.6,75) 

2.08187 

2.7728 

Optimum 
on  Plane 

(61.75,75) 

2.15123 

2.70985 

Table  1.1.  Sphere  vs.  plane,  a Steiner  problem. 
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large  region  problems  on  a sphere  via  Euclidean  geometry, 
there  is  a need  to  develop  solution  techniques  for  such 
problems . 

It  is  recognized  that  the  earth  is  a spheroid 
rather  than  a sphere.  Computation  for  the  spheroid  is 
tedious  at  best,  but  the  error  due  to  a spherical  assump- 
tion is  not  too  significant.  Use  of  spheroidal  equations 
would  only  be  warranted  where  extreme  accuracy  is  sought, 
and  when  all  other  aspects  of  the  problem  are  handled  with 
utmost  care.  If  the  earth  were  represented  by  a spheroid 
with  an  equatorial  diameter  of  25  feet,  the  polar  diameter 
would  be  approximately  24  feet  11  inches  (Deetz  and  Adams 
1945).  Maling  (1973)  cites  a study  in  which  the  average 
distance  error  when  using  a spherical  assumption  is  less 
than  1 . 07o  between  200  points  in  the  United  States. 

1 . 2 Location  Problems  Considered 

This  study  deals  with  optimizing  the  point  loca- 
tion of  a single  servicing  facility  in  a three  dimensional, 
constrained,  continuous  solution  space  for  a finite  system 
of  known  fixed  demand  points  located  on  a sphere.  The 
demand  points  are  connected  with  the  servicing  facility 
by  communications/transportation  links.  The  demand  points 
cover  a large  region,  i.e.,  a region  in  which  the  planar 
assumption  introduces  considerable  error.  A single 
equality  constraint  forces  the  solution  to  be  two  dimen- 
sional in  the  three  dimensional  space  by  restricting 
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solutions  to  be  on  the  surface  of  a sphere.  That  is,  the 
communications/ transportation  links  between  the  servicing 
facility  and  the  demand  points  are  via  great  circle  arcs 
on  the  surface  of  the  sphere  which  contains  the  demand 
points.  Costs  due  to  transportation  are  a linear  function 
of  distance.  Except  for  specific  examples,  a unit  sphere 
shall  be  assumed.  This  is  done  without  loss  of  generality 
since  geodesic  distance,  s,  is  directly  related  to  the 
radius,  r,  via  s = ro. 

Two  different  objectives  will  be  considered: 


i.  Minisum . Minimize  the  total  sum  of  costs. 
Having  M existing  facilities  located  at  known  distinct 
points  P-^,  ? 2 , •••»  PM;  a new  facility  is  to  be  located 

at  point  X.  Costs  of  a transportation  nature  are  incurred 
that  are  directly  proportional  to  the  distance  metric 
between  the  new  facility  and  existing  demand  point. 


Definition  1 . 1 Given  any  two  points  X = (x^ , *2  , • • • *xn)  and 

Y = (yj,  ,y? ....  ,yn)  in  Euclidean  En  space,  and  p > 1,  the 

£ metric  between  X and  Y is: 

P 


(Note  that  p = 1 
Euclidean  norms, 
If  W.  is 

l 

number  of  trips 


X-Y 


I V I I P 1 1/p 

i = [.\|xi-yii  ] 

p r=  1 


and  p = 2 represent  the  rectilinear  and 
respectively . ) 

the  product  of  cost  per  unit  distance  and 
per  time  period  between  X and  P.,  the  tota 
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cost  per  period  is  given  by: 


f(x)  = >:  wi|x-Pi|4 

i=l  P 


The  single  facility  location  problem  using  the  minisum 


objective  is  to  determine  the  location  of  the  servicing 


facility,  say  X*,  that  minimizes  f(X),  the  period's  total 


transportation  cost. 


For  the  large  region  location  problem  one  can 


define  p(X,P^)  as  the  shortest  great  circle  distance 


between  demand  point  and  servicing  facility  X.  This 
measurement  of  distance  is  shown  to  be  a metric  by  Blumenthal 


(1961).  The  problem  then  becomes: 


Minimize  Z W.p(P. ,X) 

X i=l  1 1 


Subject  to  | X | „ =1,  Xt:E‘ 


lere  P^cE^  are  points  on  the  unit  sphere  for  i = 1 M. 


ii.  Minimax.  Minimize  the  maximum  distance. 


That  is: 


Minimize  Z 
X l 

t ** 


Subject  to  IX-PjJ^  £ 

P 


Z i = 1 M 


where  Z is  (geometrically)  the  radius  of  a sphere  centered 


at  X,  and  P^  and  X are  as  above. 


Suppose  a known  finite  number  of  points  on  the  surface  of 


a sphere  is  given  and  it  is  desired  to  locate  a single 
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point  on  the  sphere's  surface  so  that  the  maximum  great 
circle  distance  between  this  point  and  the  given  point  is 
minimized.  The  problem  is  then  formulated  as: 

Minimize  Z 
X 

Subject  to  p(X,P.)  <.  Z 


where  X and  P.  are  as  above. 

1 

1 . 3 Applicat ion  of  Research 

This  research  has  application  in  any  large  region, 
long  range  single  facility  minimax  or  minisum  location 
problem  and,  as  will  be  pointed  out  in  Chapters  III  and  IV, 
may  be  useful  in  development  of  techniques  to  handle  the 
counterpart  mul ti f aci lity  problem.  Such  topics  as  detection 
station  placement,  naval  deployment,  location  of  interna- 
tional headquarters  or  distribution/marketing  centers,  and 
location  of  long  range  weapons  systems  fall  within  the 
purview  of  this  research.  It  would  seem  especially  perti- 
nent in  the  area  of  long  range  communications.  The  impor- 
tance in  radio  engineering  stems  from  the  fact  that  radio 
transmissions  follow  a great  circle  track. 

1 . 4 Scope  and  Limitations 

The  research  is  limited  to  consideration  of  the 


single  facility  location  problem  using  a great  circle  metric. 
The  servicing  facility  is  restricted  to  the  spherical 
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surface  upon  which  a finite  number  of  "weighted"  demand 
points  are  located.  As  cited  previously,  except  for  the 
restriction  on  location  of  the  servicing  facility,  the 
problem  is  unconstrained.  Two  criteria,  minimization  of 
total  costs  and  minimization  of  maximum  costs,  are  investi- 
gated . 

Extensive  use  of  synthetic,  projective  and  analytic 
geometry  is  made  throughout  the  research  in  development  of 
solution  procedures  and  establishing  the  theoretical  results. 

There  is  a definite  gap  in  the  theory  of  location 
with  regard  to  the  great  circle  metric.  This  research  is 
motivated  by  recognition  that  there  are  situations  in  which 
it  is  necessary  to  use  special  techniques,  notably  when 
the  regions  considered  are  larger  than  a hemisphere.  The 
research  investigates  the  properties  of  such  problems  and 
develops  approaches  for  handling  them.  It  is  hoped  this 
work  will  ultimately  stimulate  interest  in  the  development 
of  efficient  approaches  for  the  multifacility  location 
problem  on  the  sphere. 

I . 5 Order  of  Presentation 

Due  to  the  magnitude  of  work  that  has  been  done  on 
the  generalized  Weber  problem  and  other  location  models , 
Chapter  II  contains  a survey  of  only  the  literature 
directly  related  to  the  research  effort.  Extensive  bib- 
liographies on  the  general  literature  are  referenced. 

Chapter  III  presents  the  research  findings  concerning  the 
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minisum  single  facility  problem  using  a great  circle  metric. 

Two  heuristic  algorithms  are  developed  for  solution  of  the 

problem.  In  Chapter  IV  the  related  minimax  problem  is 

examined.  It  is  shown  that  existing  algorithms  for  the 
2 3 

problem  in  E and  E can  be  adapted  to  solve  the  large 
region  problem  in  certain  cases.  Chapter  V presents  compu- 
tational experience  with  the  algorithms  developed  in  Chapter 
III.  A number  of  example  problems  are  solved  and  results 
discussed.  The  research  is  summarized  and  recommendations 
for  future  research  are  made  in  Chapter  VI. 


A 


\ 

i 


This  chapter  presents  a review  of  previous  research 
pertinent  only  to  location  problems  in  large  regions.  As 
mentioned  in  Chapter  I,  two  specific  types  of  objective 
functions  are  of  interest: 

1.  Minimization  of  total  costs  (Minisum),  and 

2.  Minimization  of  maximum  distance  (Minimax). 
Concerning  the  general  literature  on  location,  one  may  refer 
to  extensive  bibliographies  (Francis  and  Goldstein  1974; 

Lea  1973) , and  a textbook  (Francis  and  White  1974) . 

Only  recently  have  researchers  shown  any  interest 
concerning  location  problems  of  the  minisum  variety  for 
regions  so  large  that  a Euclidean  (planar)  assumption  is 
not  appropriate.  Wendell  (1971)  provided  a brief  discussion 
of  the  minisum  single  facility  location  problem  on  the 
earth's  surface.  An  approximation  technique  was  given 
which  used  Schwartz's  inequality  to  find  an  explicit 
approximate  solution  (see  Section  III. 8).  Methods  for 
obtaining  an  exact  solution  were  not  considered. 

Lea  (1973)  states  that  some  work  has  been  done  on 
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the  problem  by  Krolikowski.  The  Newton-Raphson  iteration 
method  was  evidently  used  to  solve  the  problem  and  contouring 
was  proposed  both  to  supplement  problem  solution  and  to 
overcome  geographic  infeasibilities. 

Katz  and  Cooper  (1975,  1976)  discussed  the  theoreti- 
cal and  computational  aspects  of  the  problem,  with  some 
attention  paid  to  other  metrics  on  the  sphere.  An  interac- 
tive method  employing  a normalized  gradient  and  an  accelera- 
tion scheme  due  to  Steffensen  (Henrici  1964)  is  used  for 
finding  a local  optimum.  Steffensen's  technique  is  a 
standard  method  for  accelerating  convergence  in  an  itera- 
tive algorithm. 

Nothing  in  the  literature  indicates  a complete  or 
explicit  work  related  to  the  subsequent  research  effort. 
Nearly  all  the  published  work  that  refers  to  location  in  a 
large  region  or  on  a spherical  surface  has  been  concerned 
with  the  minimax  criterion.  However,  most  of  the  effort 
in  the  past  has  been  directed  toward  solving  the  problem 
for  the  continuous  case,  i.e.,  locating  n facilities  on 
the  sphere  so  that  the  maximum  distance  to  any  point  on 
the  sphere  is  minimized. 

The  attempts  at  solving  this  problem,  and  a concise 
statement  of  it,  are  given  by  Toth  (1973).  An  excellent 


review  of  solutions  for  the  cases  n = 5 and  n = 7 is  found 
in  Meschkowski  (1966).  There  are  only  a few  values  of 
n for  which  the  solution  is  known  and  for  arbitrary  values 
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of  n there  is  not  even  a reasonable  conjecture  concerning 
the  solution. 

Concerning  the  tninimax  criterion  when  only  a 

finite  number  of  demand  points  exist  on  the  sphere  (the 

problem  of  interest),  evidently  nothing  has  been  done  in 

2 

over  a century.  The  problem  in  E was  first  proposed  by 

Sylvester  (1857) . As  cited  by  Sylvester  (1860)  , the  problem 
2 

was  solved  for  E by  Peirce.  The  technique  was  redis- 
covered by  Chrystal  (1885)  . A more  recent  presentation  of 
the  approach  is  given  by  Rademacher  and  Toeplitz  (1957) 
in  a book  of  mathematical  diversions.  In  Sylvester's 
paper  the  problem  was  erroneously  claimed  to  be  completely 
analogous  to  the  one  on  the  sphere.  This  claim  is  discussed 
further  in  Section  IV. 4. 

2 

Peirce's  rudimentary  technique  for  E is  suitable 
for  solving  the  problem  by  hand.  Recently  Elzinga  and 
Hearn  (1972a)  presented  a similar  algorithm  which  was 
conducive  to  programming  on  a computer  and  hence  is  more 
efficient.  No  such  efficient  technique  exists  for  the 
counterpart  problem  on  the  sphere. 

Since  little  has  been  done  concerning  the  large 
region  location  problems,  a Euclidean  assumption  was  made, 
using  latitude  and  longitude  as  Cartesian  coordinates. 

Any  large  region  "real  world"  example  problem  in  the 
literature  is  invariably  concerned  with  not-so-large  a 
region  situated  in  the  middle  latitudes,  thus  keeping  the 


A 
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error  due  to  a Euclidean  assumption  within  bounds  of 
acceptability.  Notable  examples  of  such  problems  are 
presented  by  Kuhn  and  Kuenne  (1962),  Smallwood  (1965), 
and  Chapelle  (1969).  Two  of  these  problems  will  be  dis- 
cussed in  more  detail  in  Chapter  V. 

Kuenne  and  Soland  (1971,  1972),  in  a study  of  the 
multifacility  generalized  Weber  problem,  convert  latitude 
and  longitude  to  coordinates  on  a Mercator  Projection 
(Warntz  and  Wolff  1971),  compute  all  distances  as  rhumb- 
line map  distance,  and  convert  to  approximate  great  circle 
distances.  Although  it  worked  well  for  their  purposes,  the 
error  in  this  method  increases  dramatically  as  either  polar 
region  is  approached  or  regions  larger  than  a hemisphere 
are  considered.  Their  approach  is  discussed  further  in 
Section  III. 14.2. 


Considering  the  foregoing,  it  is  evident  that 
researchers  were  aware  of  the  problems  encountered  in 
using  a Euclidean  assumption  when  demand  points  are  scat- 
tered over  a large  region,  but  to  date  little  has  been 
done  to  resolve  them. 


CHAPTER  III 

MINISUM  SINGLE  FACILITY  LOCATION  PROBLEMS 
III.l  Introduction 

In  this  chapter  minimization  of  the  total  costs  for 
a single  facility  location  problem  on  a sphere  using  the 
great  circle  metric  is  addressed.  Two  new  algorithms  for 
solution  of  the  ''minisum"  problem  are  developed,  and  a 
number  of  properties  of  the  problem  are  presented. 

As  stated  in  Chapter  I,  a general  form  of  the  prob- 
lem for  l metrics  is  given  by: 

M 

min  f(x)  = I W. |X-P. | 0 
X i=l  1 1 p 

where:  P^  for  i=l M are  fixed  demand  points 

VL  for  i=l,...,M  are  non-negative  weights,  and 

X is  the  unknown  location  of  the  servicing  facility, 

2 

XeE^ 

The  objective,  then,  is  to  determine  the  location  of  the 
new  facility,  say  X , that  minimizes  f(X),  the  total 
"transportation”  cost.  Cost  is  considered  to  be  strictly 
a function  of  great  circle  distance  and  demand  point  weights. 
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I I I . 2 Problem  Formulation 

On  the  sphere,  formulation  of  the  problem  requires 
use  of  inverse  trigonometric  functions.  There  are  a number 
of  formulations,  each  having  its  own  advantages  and  disad- 
vantages. In  general,  the  problem  can  be  formulated  as 


M 

min  I W.p(P. ,X)  (3.2.1) 

X i=l  1 1 

Subject  to  |XL  =1,  XcE^ 

2 

where:  p(P^,X)  is  the  great  circle  distance  between  demand 

point  P.  and  servicing  facility  X,  and  | • | „ is  the 

j-  2 

Euclidean  norm.  Without  loss  of  generality,  the  problem 
is  formulated  for  the  unit  sphere  since  arc  length  on  a 
spherical  surface  is  directly  proportional  to  the  radius. 

Recognizing  that  the  shortest  distance  between  two 
points  on  a sphere  is  via  the  shorter  great  circle  arc 
connecting  them  (Lyusternik  1964)  it  can  be  shown  that 


p(P.,X)  = 2 Arcsin  — i 2 (3.2.2) 

1 2 


Thus , one  has : 


M Ip.-xL 

min  L 2W . Arcsin  ^ — 2 (3.2.3) 

X i=l  1 1 

Subject  to  | X | ^ = 1,  XcE^ 

where:  W. , P.,  and  X are  defined  as  in  Section  III.l. 
l i’ 

It  is  easily  seen  that  the  argument  for  Arcsin  is 


restricted  between  0 and  1,  and  that  the  function  is  convex 
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in  this  region.  It  is  well-known  that  a convex  function 
operating  on  a convex  function  is  likewise  convex.  Recog- 
nizing that  the  argument  is  convex  and  that  the  objective 
function  is  a non-negative  sum  of  convex  functions  it  fol- 
lows that  it,  too,  is  convex.  Unfortunately,  though,  the 
solution  space  is  not  convex,  resulting  in  a non-convex 
programming  problem.  This  is  demonstrated  graphically  in 
Figure  3.L  Although  this  form  of  the  problem  is  revealing 
with  regard  to  properties  of  the  unconstrained  problem,  it 
was  found  to  be  not  as  efficient  computationally  in  repeated 
calculations  of  objective  function  values  as  the  next 
formulation  (about  20%  slower  on  an  IBM  370/158J) . 

A basic  result  of  elementary  calculus  is  that  if 
two  lines  and  have  direction  cosines  (A^,p^,v^)  and 
(A2,p2>v2)«  and  if  0 is  the  angle  between  L-^  and  then 
cos  0 = ^^2  + yiu2  + vlv2'  Defining  ns  the  angle 
between  the  rays  in  determined  by  OP^  and  0)?,  where  0 
is  the  center  of  the  sphere,  and  observing  that  minimizing 
the  weighted  sum  of  the  arcs  (great  circle  distances)  is 
equivalent  to  minimizing  the  weighted  sum  of  the  angles 
which  subtend  the  arcs,  the  objective  is 

M 

Min  E W.ip.  (3.2.4) 

X i=l  11 

Sub j ect  to  | X | „ =1 

2 

On  a sphere  of  radius  1 with  center  at  the  origin 
3 

in  E , demand  points  = (a^,b^,c^),  i=l M and 
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X = (x^,X2,x^)  it  follows  that: 


tj; . = Arccos  (a.x,  + b.x„  + c.x0)  (3.2.5) 

1 l 1 i2  i 3 


where:  |xL  = 1.  Converting  to  spherical  coordinates 

v'2 

one  gets  the  following  formulation: 


Min  £ W.  Arccos  (a.sinC’cosA  + 
i=l  1 x 

b^sin^sinA  + c^cos$) 


Subject  to  -it  < A £ it 


(3.2.6) 


0 <.  <t>  <.  n 


where:  x-^  = sinif>cosA  >.  0 i=l M 

X2  = sin<J<sinA 
x^  = cos<J> 

This  formulation,  equivalent  to  (3.2.3),  was  found 
to  be  efficient  for  the  two  algorithms  developed  in  this 
chapter.  Using  formulation  (3.2.6),  though,  an  essenti- 
ally unconstrained  problem  with  two  bounded  variables  is 
solved.  The  bounds  are  necessary  only  because  of  periodi- 
city of  the  objective  function.  A conspicuous  disadvantage, 
which  could  affect  its  use  in  some  algorithms,  is  that  one 
knows  little  about  the  properties  of  the  objective  function. 
The  Arccos  function  is  neither  convex  nor  concave  over  the 


domain . 


Any  point  on  the  sphere  can  be  identified  by  a two- 


tuple ( , A ) where  <p  is  the  colatitude  (0  £ A <.  u)  and  A is 
the  longitude  or  meridian  (-n  < A <.  n ) . Letting  = (<f>.,A^) 
and  X = ( 0 , A ) , via  the  Law  of  Cosines  for  Elliptic  Geometry 


■ 
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(Kay  1969)  and  properties  of  the  cosine  function  it  follows 
that : 

cos[p(P^,X)]  = cos!^cos!  + s in!  ^s  in!cos  ( A - A 
and  thus: 

P(Pi,X)  = Arccos  (A^) 

The  goal,  then,  is  to: 

M 

minimize  E W.  Arccos  (A.)  (3.2.7) 

i=l  1 1 

Subject  to  0 <.  ! <.  it 

-7T  < A <.  TT 

This  formulation,  although  similar  to  the  previous 
one,  was  found  to  be  computationally  inefficient.  This  is 
primarily  due  to  time  required  for  Taylor  Series  approxi- 
mations of  trigonometric  functions  on  a digital  computer. 

Its  principal  advantages  are  that  one  works  directly  with 
spherical  coordinates,  along  with  the  property  that  the 
formulation  is  essentially  an  unconstrained  problem. 

I I I . 3 Fundamental  Properties 

In  this  section  the  non-convexity  of  the  problem 
is  demonstrated,  along  with  non-differentiability  of  the 
objective  function  and  the  fact  that  the  domain  of  the 
objective  function  is  restricted. 

I I I . 3 . 1 Non-convexity 

The  unfortunate  characteristic  of  non-convexity, 
which  is  not  a factor  in  the  planar  single  facility  location 
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problem,  is  best  seen  through  an  example.  Consider  the 
three  equally  weighted  demand  points  in  Data  Set  D1  of 
Appendix  B and  suppose  the  objective  is  to  minimize  the 
total  sum  of  great  circle  distances.  A three  dimensional 
graph  of  the  objective  function  values  as  a function  of 
and  A is  seen  in  Figure  3.1.  It  is  depicted  as  a flattened 
sphere  with  infinite  distortion  along  the  <J>A  plane  at  the 
North  and  South  Poles.  It  is  clear  from  the  graph  that 
local  minima  exist  at  P^,  > an<^  ^3>  with  possible  doubts 

as  to  what  is  happening  behind  the  peak  at  P^'.  It  turns 
out  that  all  the  local  minima  are  visible,  and  it  will  be 
seen  in  Chapter  V that  the  global  minimum  occurs  at  P^. 

Figure  3.2,  a contour  plot  of  the  minisum  objective 
function  using  Data  Set  D2  of  Appendix  B,  exhibits  the 
possible  existence  of  alternate  optima  on  the  sphere, 
namely  the  entire  arc  ^2^3'  The  8raPh  t*ie  objective 
function  over  the  great  circle  arc  as  a function  of  longitude 
is  piecewise  linear. 

I I I . 3 . 2 Non-dif ferentiability 

Besides  the  non-convexity  of  the  problem,  difficul- 
ties arise  as  to  differentiability  of  the  objective  function, 
necessitating  consideration  of  techniques  to  circumvent 
these  difficulties  during  any  search  procedure.  Although 
non-differentiability  occurs  only  at  the  demand  points 
for  the  counterpart  planar  Euclidean  norm  problem,  on  the 
sphere  the  objective  function  is  non-dif ferentiable  at  both 
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the  demand  points  and  at  the  unique  points  antipodal 
(p(x,y)  = it)  to  each  of  the  demand  points.  Referring  to 
Figure  3.1  again,  the  "knife-edged  surface"  described  by 
Vergin  and  Rogers  (1967)  appears  at  these  points.  The 
reason  for  this  occurrence  can  be  seen  by  examining  the 
derivative  of  Arccos  X: 

SArccos  X , o 

9X  Wl-X^ 

Note  that  the  first  derivative  is  undefined  when  X = ^ 1 , 
which  occurs  for  arcs  of  length  0 and  tt  , which  in  turn 
correspond  to  the  situations  when  one  is  at  a demand  point 
or  its  antipodal  point. 

Concerning  the  search  techniques  which  will  be 
developed,  one  of  them  requires  consideration  of  this 
property  since  it  employs  a gradient  search.  The  other 
algorithm  is  derivative-free. 

f 

IV. 3. 3 Restricted  Domain  of  Objective  Function 

Another  property  of  the  single  facility  minisum 
sphere  problem  not  encountered  in  the  planar  case  is  that 
the  domain  for  the  objective  function  is  restricted.  This 
is  because  the  real  argument  Z of  the  Arccos  function  may 
not  exceed  the  bound  |z|  <.  1.  Given  that  X is  a possible 
servicing  facility  location,  certainly  the  objective  func- 
tion is  defined  for  all  possible  X contained  within  or  on 
the  unit  sphere,  but  not  defined  for  all  X outside  the 


sphere.  In  Figure  3.2  this  is  demonstrated  for  a four  point 
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problem  in  the  degenerate  case  of  location  on  a great 
circle  arc.  The  "outer"  area  is  not  within  the  domain  of 
definition.  As  the  number  of  demand  points  increase  around 
the  sphere,  note  that  the  domain  of  definition  will  shrink 
to  approximate  the  sphere. 

This  restricted  domain  of  definition  must  be  con- 
sidered when  using  any  search  technique  calling  for  projec- 
tions. Care  must  be  taken  to  insure  that  the  search  does 
not  leave  the  domain  of  definition,  thus  causing  premature 
termination . 

I I I . 4 Dominance  Properties  on  the  Sphere 
I I I . 4 . 1 Introduction 

Although  the  general  minisum  problem  on  the  sphere 
has  many  undesirable  properties,  there  are  situations  in 
which  one  can  reduce  the  search  region  for  a global  optimum. 
For  example,  it  seems  plausible  and  intuitive  that  if  all 
the  demand  points  are  located  on  a hemisphere,  then  any 
search  for  an  optimal  solution  can  be  restricted  to  this 
region.  Using  specific  results  of  convexity  theory  for 
spherical  geometry,  this  intuitive  concept  is  now  generalized 
to  demonstrate  that  any  search  for  an  optimal  solution  to 
(3.2.1)  can  be  restricted  to  the  spherically  convex  hull  of 
the  demand  points. 

The  major  results  of  this  section  are  based  upon  a 
generalization  of  Kuhn's  (1967)  concept  of  dominance  due 
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to  Wendell  and  Hurter  (1973).  It  is  described  by: 

Definition  3.4,1  A point  x'  dominates  a point  x with 

respect  to  P-^ , . . . , P^ P^ . demand  points,  if  and  only  if 

p(x',Pi)  £ pCx.P^  for  all  i. 

As  an  immediate  result,  if  x'  dominates  x with 

respect  to  P^  for  all  i,  and  is  a non-negative  constant 

for  all  i,  x'  has  the  property  that: 

M M 

I W . p (x ' , P . ) < l W.p(x,P.) 
i=l  1 1 ~ i=l  1 1 

So , by  showing  for  any  x^V,  where  V is  a spherically 
convex  hull,  that  there  exists  an  x'eV  such  that  x'  domi- 
nates x one  is  assured  of  the  existence  of  an  optimal  solu- 
tion x*  such  that  x*  is  in  V. 

Prior  to  establishing  the  main  results  it  is  neces- 
sary to  introduce  some  concepts  and  lemmas.  Terminology  and 
proofs  of  supportive  lemmas  are  in  Appendix  A. 

I I I . 4 . 2 Dominance  within  Spherically  Convex  Hull 

In  this  section  it  is  established  that  any  search 
for  an  optimal  solution  to  (3.2.1)  can  be  restricted  to  the 
spherically  convex  hull  containing  the  demand  points. 

Theorem  3.4.1  Given  a set  of  demand  points 
{PjJi=l M}  (not  consisting  of  two  antipodal  points) 

2 

whose  convex  hull  is  V = convtP^ | i=l , . . . ,M} . For  any  xeS 

2 

such  that  x^V,  there  exists  x*cVcS  such  that  x*  dominates 
x in  the  great  circle  metric  with  respect  to  P^  , i=l M. 
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Proof : By  Lemma  A. 8,  V is  closed.  So  applying 

Lemma  A. 4,  consider  the  support  line  L orthogonal  to  line 
xp  at  p.  The  proof  is  in  two  parts. 

Part  A.  Suppose  that  V lies  entirely  on  a line 
through  x,  and  further  that  the  measure  of  V = pq  is  £ tt 
(note  that  if  the  measure  is  > tt  , then  it  is  equal  to  2tt). 

For  x^V,  designate  p as  the  closest  point  of  V to  x. 

If  p(x,p)  .>  tt  / 2 , choose  x*  as  the  midpoint  of  pq\  Since 
p(p,q)  £ tt  , then  p(x*,P^)  £ tt/2  for  all  P.cV,  i=l,...,M; 
and  x*  dominates  x. 

If  p(x,p)  < tt / 2 , consider  the  following.  If 
p(p,q)  £ p (x,p)  then  the  choice  x*  = p obviously  dominates 
x since  for  any  demand  point  P^  , (xpP^q)  holds  and  p(x,v,P^)  = 
p(p,pi)  £ p (p , q)  £ p(x,p)  £ p(x,P.) . Otherwise,  take 
x*tpq^  such  that  (xpx*)  and  p(x,p)  = p(p,x*).  Also,  define 
x'  and  p’  as  the  antipodal  opposite  points  of  x and  p, 
respectively.  Let  P^  be  a demand  point. 

Case  I.  P^epx*  (Figure  3.3a) 

Then  (xpP^x*)  holds  and  p(x*,P^)  £ p(x*,p)  = p(x,p)  £ p(x,P^) 
by  properties  of  betweenness;  hence  x*  dominates  x. 

Case  II.  P . ex*x ' (Figure  3.3b) 

Then  (xpx*P^)  holds  and  p(x,P^)  = p(x,x*)  + p(x*,P^)  £ p(x*,P^) 
by  properties  of  betweenness,  and  dominance  follows. 

Case  III.  P^ex'p'  (Figure  3.3c) 

Then  (xp'P^)  holds,  as  does  (x*x ' P^p ' ) . Now  since, 

P(p.p')  = tt  = p(p,x)  + p(x.p') 
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= p (p ,x*)  + p (x*,p ' ) 

= p (p , x)  + p (x*,p ' ) , 
it  follows  that  p(x,p')  = p(p',x*). 

Adding  p(p'.P^)  to  both  sides,  and  utilizing  the  fact  that 
(xp ' ) , one  has 

p (x , P^)  = p(p',x*)  + p (p ' ,Pi)  2l  p(p'.x*) 

Since  (x*P^p'),  it  follows  from  betweenness  properties  that 
p(p',x*)  = p(x*,Pi)  + p(Pi(p')  >_  p(x*,Pi) 

Thus  p(x,P.)  _>  p(x*,Pi)  and  x*  dominates  x. 

Since  all  possible  cases  have  been  considered,  the 
proof  of  Part  A is  complete. 


(a)  (b)  (c) 


Figure  3.3.  Theorem  3.4.1,  Part  A. 


Part  B.  Next,  assume  V does  not  lie  on  a line 
through  x.  By  Lemma  A. 5,  there  exist  exactly  two  distinct 
rays  xb  ^ and  xb2  with  h^ , b,  r bd  V such  that  the  lines 
xb^  and  xb^  are  support  lines  of  V and  each  ray  xv  for 
vt V either  coincides  with  xb^ , or  xb^ , or  lies  between 
xbj  and  xB^.  In  particular,  xp  is  between  xb , xb^  or 
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coincides  with  them,  so  the  sum  of  the  angles  ^ pxb , and 
^[pxb  , is  no  greater  than  u , and  at  least  one  of  them,  conse- 
quently, can  be  no  greater  than  it/2.  We  may  assume  without 
loss  of  generality  $ pxb^  _ ti/2.  Now  either  -$.pxb7  £ tt/2 
or  pxb0  > v/2.  The  cases  are  considered  separately. 

Case  I.  Suppose  pxb,,  £ it/ 2 (see  Figure  3.4a). 

Let  x*  = p,  and  suppose  v is  any  point  of  V not  on  line  1<p 
(there  must  exist  at  least  one  such  point  under  the  assump- 
tion)  . If  v is  on  the  b^-side  of  xp  (k  = 1,2)  then  since 
xv  coincides  with  xb^  or  lies  between  xb  ^ and  xb^,^tpxv 
£ •<£.  b^xp  r.  / 2 . On  the  other  hand,  with  v on  L or  on  the 
opposite  side  of  L as  x,  one  of  the  opposing  rays  from  p 
on  L either  coincides  with  pv  or  lies  between  px  and  pv  so 
that  •£.  xpv  ^ n/2.  Hence  in  all  cases  pxv  2 ■<£  xpv . If 
equality  holds,  p(p,v)  = p(x,v).  If  inequality  holds, 
using  the  properties  that  the  largest  side  of  a spherical 
triangle  is  opposite  the  largest  angle  (Kay  1969,  Theorem 
31.5),  it  follows  that  p(p,v)  < p(x,v).  That  is,  p(x*,v)  £ 
p(x,v)  for  all  v in  V not  on  line  xp . However,  by  connected- 
ness of  V this  holds  for  all  points  of  V.  So,  for  all 
demand  points  P . it  follows  that  p(x*,P^)  £ p(x,P^).  Domi- 
nance of  x*t V has  thereby  been  proved  for  this  case. 

Case  II.  Suppose  pxb.,  N ?r/2.  According  to  Lemma 
A. 6 one  may  choose  ut  xb.,  and  x'ebd  V such  that  p(x,u)  = 


p(u.x')  and  ux'  is  a line  of  support  of  V.  (Figure  3.4b). 

In  this  case  put  x*  = x'.  Let  u'  be  the  antipodal  point  of  u. 
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Then  both  triangles  xux'  and  xu'x'  are  isosceles.  Once 
again  xx'  either  coincides  with  xb^  or  falls  between  xb^  and 

- — 'y  ^ ^ 

xb0.  Let  v be  any  point  of  V not  on  xx' ; then  xv  coincides 
with  xb^,  k = 1,2,  or  lies  between  them.  In  the  case  when 

xv  lies  on  xb*,  or  between  xx ' and  xb^ , ray  x'v  lies  between 

. — ^ r 

x'u  and  x'y  (y  is  any  point  on  xx ' such  that  (xx'y))  or 
coincides  with  x'u.  Hence  -^xx'v  > ^ xx'u  =4.  x'xu  £ 

^ vxx1  and  by  properties  of  triangles  having  unequal  angles, 
p(x',v)  £ p(x,v).  The  same  result  holds  when  lies  on 
the  other  side  of  line  xx ' , appealing  to  equal  angles 


4u'xx'  an 


d u'x'x.  Again,  as  in  Case  I,  p(x*,v)  £ p(x,v) 


holds  for  all  veV  and  specifically  for  all  demand  points 
= veV,  proving  dominance  of  x*  in  this  last  case. 
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Corollary  3.4.1  For  the  great  circle  metric,  there  exists 
an  optimal  solution  x*  to  the  problem  (3.2.1)  such  that 
x*  e conv  { P . | i=l , . . . ,M}  = V. 

Proof : By  using  a spherical  coordinate  system  the  con- 

straint is  satisfied  implicitly.  The  result  follows  from 
the  observation  that  if  Xq  dominates  x,  then 

M M 

E WiP(x  P ) < E Wp(x,P  ). 
i=l  i=l 

One  need  only  consider  elements  within  the  region  since  if 
there  exists  Xq  outside  of  V which  is  optimal,  Theorem 
3.4.1  guarantees  the  existence  of  an  x*  within  V which 
dominates  Xq . 

I I I . 4 . 3 Determining  Demand  Points'  Hull  Characteristics 
In  order  to  apply  the  foregoing  results , it  is 
of  course  necessary  to  insure  that  the  set  of  demand  points 
can  be  contained  in  a hemisphere.  For  if  this  is  not  possi- 
ble, the  spherically  convex  hull  is  the  entire  sphere,  and 
the  search  region  is  not  reduced  at  all.  The  logical 
procedure  is  to  coordinatize , plot  the  points  on  a sphere 
and  visually  verify  whether  the  set  is  containable  in  a 
hemisphere . 

Blumenthal  (1956)  defined  a global  subset  as  a 
2 

subset  G of  S which  is  not  contained  in  any  hemisphere 
2 

of  S . For  the  case  of  m=4  points,  he  established  neces- 
sary and  sufficient  conditions  for  the  points  to  form  a 
global  subset.  Clearly  for  m £ 3,  the  points  can  be  contained 
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in  a hemisphere.  Nothing  is  known  about  the  general  case 
of  m > 4 . 

What  follows  is  a rudimentary  heuristic  procedure 
for  determining  whether  the  set  { P . | i=l , . . . ,m}  forms  a 
global  subset.  First,  rank  all  inter-demand  point  distances 
from  the  smallest  to  the  largest  and  select  the  most  distant 
pair.  Using  either  of  the  points,  say  , transfer  the 
pole  to  that  point.  Rank  the  longitudes  of  all  the  remain- 
ing demand  points  in  ascending  order,  and  see  if  any  abso- 
lute difference  between  the  longitudes  of  "adjacent"  demand 
points  (with  reference  to  rank  order  of  longitudes)  is 
greater  than  or  equal  to  u . If  so,  one  of  the  two  hemi- 
spheres determined  by  the  line  through  P^  and  the  point 
having  minimum  of  the  two  longitudes  contains  the  set  of 
demand  points.  It  is  the  hemisphere  which  contains  the 
point  with  the  maximum  of  the  two  longitudes.  The  absolute 
difference  between  the  minimum  and  maximum  longitude  of  all 
demand  points  (other  than  P^)  is  given  by  2rc  - (max  longitude 
- min  longitude) . 

If  no  absolute  difference  is  greater  than  ti  , 
choose  the  other  point  of  the  original  pair,  unless  it  is 
an  antipodal  point,  and  repeat  the  procedure.  This  is 
continued  through  the  list  of  inter-demand  point  distances 
until  either  an  enclosing  hemisphere  is  found  or  all  points 
have  been  considered. 


it. 
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III.  A. 4 Summary 

In  this  section  it  was  shown  that  the  search  for  an 
optimal  solution  can  be  restricted  to  the  convex  hull  of 
the  demand  points . Note  that  if  the  set  of  demand  points 
is  not  containable  in  a hemisphere,  the  convex  hull  is  the 
entire  sphere,  and  no  reduction  in  the  search  area  is  possi- 
ble. The  only  sets  of  demand  points  of  interest,  then,  are 
those  containable  in  a hemisphere. 

A rudimentary  procedure  is  provided  to  determine 
whether  the  spherically  convex  hull  is  containable  in  a 
hemisphere  for  a particular  set  of  demand  points.  Cer- 
tainly, based  upon  these  results,  one  would  not  pick 
starting  solutions  outside  the  convex  hull  when  applying 
an  iterative  technique. 

IV .  5 Analogue  Models 

Mechanical  and  electrical  analogue  models  for  the 

Euclidean  norm  single  facility  location  problem  have  been 

devised  and  successfully  employed.  The  purpose  of  this 

section  is  to  extend  these  models  to  handle  the  related 
2 

problem  on  S . 


» 


' 

I 

I 


I I I . 5 . 1 Mechanical  Analogue 

A description  of  the  basic  analogue  model  is  given 
by  Eilon  et  al.  (1971),  Francis  and  White  (1974),  Lyusternik 
(1964),  and  Haley  (1962),  among  others.  Utilizing  strings 
and  weights,  it  was  introduced  by  Georg  Pick  in  the  early 
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1900's  in  Alfred  Weber’s  Uber  den  Standort  der  Industrien 
(1909).  Description  of  a case  study  using  this  model  is 
given  by  Burstall  et  al.  (1962). 

The  model  which  follows  will  solve  the  single 
facility  great  circle  metric  problem  (3.2.1).  It  can  be 
used  for  regions  as  large  as  a hemisphere  without  diffi- 
culty . 

First,  a highly  polished  sphere  (to  minimize  fric- 
tion effect)  is  coordinatized  and  demand  points  are  plotted. 
Holes  are  then  drilled  at  the  demand  point  locations  and 
strings  are  passed  through  the  holes  with  the  ends  tied  to 
a small  ring  resting  on  the  exterior  surface  of  the  sphere. 
The  other  ends  are  passed  through  a stationary  ring  at 
the  sphere's  center.  Weights  proportional  to  the  respec- 
tive demands  at  each  point  are  tied  to  the  appropriate 
strings  (see  Figure  3.5). 


Figure  3.5.  A mechanical  model. 
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The  outer  ring  (on  frictionless  bearings)  is  pulled 
to  one  side  and  then  released.  The  weights  will  pull  the 
ring  to  a point  of  minimum  potential  energy  at  which  there 
is  a local  minimum  objective  function  value.  Due  to  the 
demonstrated  non-convexity  of  the  problem,  it  is  possible, 
although  unlikely,  that  the  ring  could  start  at  a local 
maximum  and  not  move.  For  this  reason,  it  should  be  moved 
slightly  to  one  side  to  see  if  it  returns  to  the  same  point. 
Local  maxima  and  minima  demonstrate  points  of  unstable  and 
stable  equilibria,  respectively.  The  unstable  equilibria 
occur  due  to  the  duplicity  of  paths  between  any  two  points; 
that  is,  the  long  and  short  great  circle  arcs.  This  property 
is  treated  in  detail  by  Lyusternik  (1964) . 

In  regions  larger  than  a hemisphere  physical  diffi- 
culties can  arise  in  the  model  as  the  ring  approaches  a 
point  y antipodal  to  a demand  point  x (p(x,y)  = ti  ) . At 
such  a point  there  are  an  infinite  number  of  paths  of  equal 
length  to  the  antipodal  demand  point.  Yet,  by  moving  an 
arbitrarily  small  distance  in  any  direction,  the  shortest 
path  becomes  unique.  The  difficulty  arises  in  the  physical 
slipping  of  the  string  to  a position  180°  opposite  to  the 
existing  position  just  prior  to  the  ring's  passing  through 
the  antipodal  point.  Friction  and  interference  from  other 
strings  would  prevent  free  movement  of  the  ring  in  such  a 
case . 


The  advantage  of  this  model  is  that  rough  estimates 
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of  location  for  local  optima  can  be  found  visually  rather 
quickly  by  starting  the  ring  at  various  places  around  the 
sphere,  thus  providing  good  starting  solutions  for  any 
iterative  technique. 

Disadvantages  include  the  friction  effect,  which 

may  make  the  final  position  of  the  ring  indeterminate.  Due 

to  the  curved  surface  and  necessity  of  the  ring  at  the  center 

of  the  sphere  it  can  be  expected  that  friction  will  have  a 

2 

greater  effect  than  in  the  corresponding  model  for  E . 

Also,  the  non-convexity  of  the  problem  results  in  the  pos- 
sibility of  only  getting  local  minima  and  missing  the  global 
optimum. 

As  mentioned  by  Eilon,  et  al.  (1971),  perhaps  the 
greatest  disadvantage  is  the  fact  that  the  method  does  not 
evaluate  the  cost  function.  The  next  model  resolves  this 
difficulty  and  essentially  eliminates  the  above  mentioned 
disadvantages . 

1 1 1 . 5 . 2 Electronic  Analogue 

Hitchings  (1969)  recognized  that  most  of  the 

literature  on  location  problems  to  that  date  was  confined 

to  problems  in  the  planar  state.  He  addressed  the  solution 
3 

of  problems  in  E and  developed  an  electronic  analogue  model 

to  solve  these  problems. 

Hitchings'  model  can  effectively  be  modified  to 
2 

solve  the  problem  in  S space  using  the  great  circle  metric. 
One  must  utilize  a non-conductive  sphere,  set  up  as  in  the 
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mechanical  analogue  model.  The  string,  however,  is  replaced 
by  resistive  wire  in  which  distance  is  proportional  to  the 
length  of  the  resistor.  Weights  are  needed  solely  to  insure 
that  the  wire  conforms  to  a geodesic  path  on  the  sphere. 

The  varying  weights  of  the  demand  points  are  accommodated 
via  changes  in  resistivity  of  material  or  cross-sectional 
area  of  the  wire.  As  pointed  out  by  Hitchings,  such  a 
model  can  even  handle  nonlinear  costs  by  segmenting  the 
wires  into  appropriate  lengths  and  varying  resistances. 

Based  upon  a slide-wire  concept,  the  model's  cir- 
cuit diagram  (Hitchings  1969)  is  similar  to  a Wheatstone 
bridge  circuit.  The  difference  lies  in  the  fact  that  the 
objective  is  to  maximize  current  flow  rather  than  find  a 
null  point.  Note  that  for  a fixed  voltage,  maximization 
of  current  flow  is  equivalent  to  minimization  of  overall 
resistance,  which  in  turn  is  the  analogue  of  the  objective 
function  for  the  single  facility  location  problem  on  the 
sphere . 

Attaching  a pen  to  the  ring  on  the  sphere's  surface, 
isocost  lines  may  be  easily  plotted  by  moving  the  pen  while 
carefully  maintaining  a constant  current  flow.  In  this 
manner  contours  can  be  plotted  to  reveal  local  minima. 

From  this  information  it  is  simple  to  determine  the  location 
of  the  global  optima  (point(s)  of  maximum  current  flow) . 


The  advantages  of  this  model  over  the  previous  one 
are  obvious.  Relative  cost  function  values  are  available 
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via  current  flow  values,  friction  has  no  effect,  and  with 
a little  care  problems  for  regions  larger  than  a hemisphere 
can  be  handled.  Since  contours  are  available,  the  model 
may  be  useful  when  considering  practical  problems  with 
geographic  infeasibilities. 

I I I . 6 St einer's  Problem  and  Fagnano's  Result 
III. 6.1  Introduction 

Consider  the  problem  of  determining  for  a triangle 
in  the  plane  the  point  at  which  the  sum  of  distances  from 
the  point  to  the  vertices  is  minimized.  This  is  a special 
case  of  the  single  facility  location  problem  where  all 
weights  are  equal.  It  is  known  that  if  any  angle  of  the 
triangle  equals  or  exceeds  120°,  the  optimum  occurs  at  the 
vertex  of  the  obtuse  angle.  If  all  angles  are  less  than  or 
equal  to  120°,  the  optimum  is  interior  to  the  triangle  at  a 
point  at  which  each  side  of  the  triangle  subtends  an  angle 
of  120°.  A recent  elegant  proof  of  the  former  property 
has  been  given  by  Sokolowsky  (1976).  Proofs  of  both  proper- 
ties are  provided  by  Courant  and  Robbins  (1941),  among 
others . 

Called  Steiner's  Problem  by  many,  it  has  a long 
history  which  is  succinctly  outlined  by  Cooper  (1963). 

Also  mentioned  by  Cooper  is  a result  due  to  Fagnano  in 

1775  showing  that  the  point  for  which  the  sum  of  the  distances 

from  the  vertices  of  a quadrilateral  is  a minimum  is  given 
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by  the  intersection  of  the  diagonals.  The  purpose  of  this 
section  is  to  consider  extension  of  the  above  results  to 
the  problem  on  the  sphere. 

I I I . 6 . 2 Steiner's  Problem  on  the  Sphere 

Consider  three  demand  points  on  the  sphere  with 
equal  weights.  Suppose  an  optimal  solution,  local  or 
global,  is  known.  Using  the  optimal  point  as  pole,  trans- 
form the  points  to  the  plane  with  an  Azimuthal  Equidistant 
Projection  (Deetz  and  Adams  1948).  Such  a projection  pre- 
serves distance  and  bearing  to  other  points  from  the  polar 
point . 

The  image  of  the  polar  point  transformed  from  the 
sphere  is  globally  optimal  on  the  plane.  Suppose  it 
were  not  optimal  on  the  plane.  Since  the  problem  is 
convex  (Love  1967)  , a move  could  be  made  over  an  arbi- 
trarily small  distance  6 > 0 in  some  direction  and  yield 
an  improvement  in  the  objective  function  value.  Since  it 
can  be  shown  that  the  objective  function  value  on  the  plane 
is  greater  than  value  on  the  sphere  for  the  corresponding 
point  (see  Section  III. 10.2)  there  would  be  an  improvement 
by  moving  in  the  same  direction  and  distance  on  the  sphere. 
This  contradicts  the  fact  that  a local  or  global  optimum 
is  reached  on  the  sphere. 

Since  the  point  is  global  on  the  plane,  it  is 
known  that  all  angle  measures  subtended  from  the  triangle 
sides  are  _>  120°,  with  strict  inequality  holding  only  when 
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the  point  is  a vertex.  The  projection  used  preserves 
angles  at  the  polar  point  so  the  property  is  retained  at 
the  inverse  image--which  is  a local  or  global  optimum  on 
the  sphere. 

It  is  important  to  note  this  result  is  a necessary 
condition  for  global  optimality  on  the  sphere.  It  is  sus- 
pected, but  not  proven,  that  if  the  points  are  such  that  no 
side  of  the  spherical  triangle  is  of  length  ^ / 2 , then  if 
any  vertex  angle  is  120°  or  larger,  its  vertex  is  the  global 
optimum.  Otherwise,  the  global  optimum  is  at  the  interior 
point  at  which  each  side  of  the  triangle  subtends  an  angle 
of  120°. 

I I I . 6 . 3 Fagnano's  Result  on  the  Sphere 

Suppose  four  equally  weighted  points  determining 
a convex  quadrilateral  (Kay  1969)  are  on  an  open  hemisphere 
and  the  optimal  point  is  X *.  The  optimal  point  occurs  at 
the  intersection  of  the  small  great  circle  arcs. 

Given  the  optimal  point,  project  to  the  plane  via 

the  Azimuthal  Equidistant  Projection.  Distance  and  bearing 

are  preserved  from  point  X on  the  plane  corresponding  to 

X *. 
s 

On  the  plane  the  optimal  is  determined  by  intersec- 
tion of  the  diagonals.  The  optimal  point  must  coincide 
with  Xp . For  if  Xp  is  not  optimal  on  the  plane,  then  improve- 
ment can  be  made  by  moving  an  arbitrarily  small  distance 
6 > 0 in  some  direction  from  X to  a point  X* . Now,  as 
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stated  in  the  previous  section,  the  objective  function 
value  of  X'  for  the  Euclidean  norm  is  strictly  greater  than 
the  objective  function  value  obtained  for  the  great  circle 
metric  using  the  point  corresponding  to  X*  for  the  problem 
on  the  sphere. 

This  means  there  is  a point  on  the  sphere  dominat- 
ing the  global  optimum--a  contradiction. 

So  Xp  and  Xg*  are  optimal  on  the  plane  and  sphere 
respectively--and  are  corresponding  points  via  the  transforma- 
tion . 


Through  the  properties  of  the  projection,  the  angles 

and  distances  are  preserved  at  X and  X *.  It  follows 

r p s 

immediately  that  the  desired  property  of  the  global  solution 


at  Xp  extends  to  the  sphere's  global  solution  as  a necessary 
condition . 


I I I . 7 A Conjecture  Concerning  Global  Optimality 

In  working  with  the  problem  on  the  sphere,  and 
after  plotting  a number  of  example  problems  three  dimen- 
sionally, an  interesting  characteristic  came  to  light. 
Although  not  established  theoretically,  it  is  important  to 
report  since  its  validity  would  significantly  impact  on 
conclusions  concerning  global  or  local  optimality. 


Conj ecture  3.7.1  If  (1)  all  demand  points  for  problem 
(3.2.1)  are  located  within  an  octant  of  the  unit  sphere,  or 
a disk  of  diameter  < tt/2  and  (2)  at  least  three  of  the 
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demand  points  are  non-collinear , then  the  objective  func- 
tion is  unimodal  within  the  region. 

The  conjecture  incorporates  two  fundamental  assump- 
tions. First,  that  there  exists  a unique  global  minimum 
within  the  region.  This  is  assured  by  requiring  non- 
collinearity  of  demand  points.  Second,  the  characteristic 
of  the  objective  function  under  this  assumption  is  unimodal. 
An  explicit  discussion  of  generalized  unimodality  in  n- 
dimensions  is  provided  by  Sivazlian  and  Stanfel  (1975)  . A 
function  f is  said  to  be  un imodal  over  a region  S if  there 
exists  a path  from  xeS  to  the  global  optimum  x*  over  which 
f is  strictly  decreasing. 

If  this  conjecture  is  true,  then  a convergent 
algorithm  searching  over  the  region  will  result  in  a 
global  optimum.  This  follows  since  the  dominance  results 
of  Section  III. 4. 2 permit  the  search  for  a global  to  be 
restricted  to  the  spherically  convex  hull  of  the  demand 
points.  The  conjecture  eliminates  the  possibility  of  having 
local  minima  within  the  region. 

A typical  contour  of  a problem  with  the  given 
properties  is  found  in  Figure  3.6.  It  plots  the  objective 
function  values  for  the  six  point  unequal  weight  minisum 
problem  using  the  Data  Set  D3  of  Appendix  B. 

Concerning  the  disk  of  diameter  tt  ^ 2 and  the  octant, 
both  regions  can  contain  sets  of  demand  points  not  contained 
in  the  other.  The  octant,  also  known  as  the  Reuleaux 


a.  e 
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Figure  3.6.  3-D  minisum  functional  plot;  Data  Set 

D3,  Appendix  B. 
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Triangle  on  the  sphere,  has  the  smallest  area  of  all  figures 
of  equal  altitude.  Using  techniques  of  Chapter  IV  it  is 
possible  to  determine  whether  the  demand  points  are  located 
in  a disk  of  diameter  <_  it / 2 . It  is  not  so  easy  to  determine 
analytically  whether  all  points  are  contained  with  an 
octant.  However,  one  can  always  plot  the  points  on  a 
sphere  and  use  an  overlay. 

It  is  suspected  that  the  unimodality  property  is 
related  to  the  fact  that  many  Euclidean- like  properties 
hold  in  such  regions  (Kay  1969) . These  properties  hold 
since  no  distance  between  any  demand  points  or  between  a 
possible  location  and  any  demand  point  is  greater  than  or 
equal  to  tt/2.  In  such  regions,  for  example,  it  is  known 
that  the  hypotenuse  of  a spherical  right  triangle  is  the 
longest  side.  This  property,  among  others,  does  not  hold 
in  the  general  case. 

I I I . 8 An  Approximate  Solution 

Wendell  (1971),  formulating  the  problem  as  (3.2.3), 

approximated  the  objective  function  by  recognizing  that 

2 

Arcsin  y & Try  /2  for  ye  I 0 , 1 1 . Through  Schwartz's  inequality, 
the  solution  to  the  revised  problem  is  seen  to  be 
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I 

Wendell  offers  x^*  as  an  approximate  solution  to  (3.2.1). 

Through  simple  algebraic  manipulation  it  can  be 

3 

shown  that  if  the  centroid  is  found  in  E and  projected  to 
the  spherical  surface  (normalized  in  the  case  of  the  unit 
sphere),  the  resulting  point  is  equivalent  to  Xq* . That 
is,  Wendell's  approximate  solution  is  the  projected  centroid. 

As  will  be  seen  in  Chapter  V,  it  is  not  difficult 
to  construct  examples  in  which  the  projected  centroid  is 
far  from  the  optimal  solution,  nor  can  one  always  success- 
fully utilize  it  as  a starting  solution  in  any  iterative 
technique  in  order  to  find  a global  optimum.  It  is  not 
without  value,  however,  as  will  be  discussed  later. 

I I I . 9 Bounds  for  Unconstrained  Objective  Function 

In  this  section  bounds  are  generated  for  the  uncon- 
strained great  circle  metric  single  facility  location  problem. 
The  results  are  an  extension  of  those  due  to  Pritsker  and 
Chare  (1970)  for  the  Euclidean  problem. 

The  bounds  are  based  upon  optimal  solutions  for  the 

rectilinear  and  Euclidean  norm  single  facility  problem  in 
3 

E . As  mentioned  earlier,  efficient  procedures  exist  for 


M i o M 99  M o 9 
V(Z  w.Pj)Z+(  Z w.Pt)Z+(  Z w.Pf)Z 
i=l  11  i=l  1 1 i=l  1 1 


r 


solving  these  problems. 

The  bounds  are : 


[R2(xR)+R2(yR)+R2(zR) ]%  < A (X°)  < A(XR)  (3.9.1) 


E(XE)  < A (X°)  < A(XE) 


(3.9.2) 


R R R R 

where:  X = (x  ,y  ,z  ) is  the  optimum  rectilinear  solution, 

E E E E 

X = (x  ,y  ,z  ) is  the  optimum  Euclidean  solution, 

X°  = (x^,y^,z^)  is  the  optimum  great  circle  metric 
solution, 

and  R(x) , E(x) , and  A(x)  are  the  objective  function  values 
for  the  rectilinear,  Euclidean  and  great  circle  metric 
problems,  respectively. 

Note  that  although  the  bounds  in  (3.9.2)  are  more  "costly" 
to  compute,  they  will  in  general  be  the  tighter  of  the  two. 
This  becomes  evident  in  the  derivations  which  follow: 

The  right  hand  inequalities  are  both  clearly  true 
for  the  optimum  of  the  objective  to  minimize  A is  X^ , so 
A (X°)  < A(X)  for  all  X. 

In  particular,  A(X^)  <_  A(XR)  and 

A (X°)  < A(XE) . 

The  fact  that  A(X^)  E(X^)  is  intuitive,  but  can  be 

established  geometrically.  Concerning  the  relationship 
^ \ 

between  arc  ANB  and  chord  ASB  (see  Figure  3.7),  note  that 
Z = \ ANB,  and  SB  = % ASB. 
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Figure  3.7.  Relationship  between  chords  and  arcs. 

Now,  without  loss  of  generality  <f>  can  be  limited  to 
the  region  [ 0 , tt  / 2 ] . Consider  the  right  triangle  NSB  where 
£ NSB  = tt/2.  Since  the  hypotenuse  of  a right  triangle  is 
the  longest  side,  NB  _>  SB.  The  shortest  distance  between 
any  two  points  is  a straight  line  so  NB  >_  NB,  and  it  fol- 
lows that  Z = NB  £ SB.  So  ASB  £ ANB . 

Since  A(X^)  is  a weighted  sum  of  great  circle 
arcs  of  form  ANB  and  E(X^)  is  the  corresponding  Euclidean 
sum,  A(X°)  > E(X°) . 

Clearly  E(XR)  £ E(X^).  Otherwise  the  optimality  of 

g 

XJ  for  the  Euclidean  norm  is  contradicted.  This  estab- 
lishes (3.9.2) . 

For  (3.9.1),  application  of  Minkowski's  inequality 
establishes  that  E(XR)  £ tR2(xR)  + R2(yR)  + R2(zR)]^  and 
the  bound  follows  via  transitivity. 

Since  X is  optimum  for  the  rectilinear  problem  one 

has  that 


[R2(xR)+R2(yR)+R2(zR) ]%  < 


[ R2(xE)+R2(yE)+R2(zE)  ]* 


A 
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Now, 

[ R2  (xlJ)+Ri‘  (y^j+R^  (z^) 


2F  ? F 2 F 3-  M 

*-  / .-L*  \ i r>  ^ \ _ 


M P o M p 

E w.  | x -a . i)  + ( £ w.  |y  -b.i)' 
i=l  1 1 i=l  1 1 


+ ( £ wt|  zK-cil  )2]^ 
i=l  1 1 


which,  by  Minkowski's  inequality  is 
< -n/(w1  | xE-ajJ  )2+(w1|  yE-b  L | j^+CwJ  zE-cp  2 + ... 


+-/(^M!xE-aMl)2  + (wMlyE-bMD2+(wMizE-cMl)"2 


V /Tl  x2,  , E , N 2 , , E 71 
= E w.  -J  (x  -a.)  +(y  -b.)  + (z  -c.) 

i=l  1 1 1 1 

= E(XE) 


„ , E E 
E (x  , y 


E 

z 


) 


III.  10  A Planar  Projection  Algorithm  (PPA) 

III. 10.1  Introduction 

In  this  section  the  first  of  two  algorithms  for  the 
minisum  single  facility  location  problem  is  developed. 

It  capitalizes  on  existing  solution  techniques  for  the 
planar  case  and  makes  use  of  the  fact  that  isometric 
(distance-preserving)  transformations  from  the  sphere  to 
the  plane  from  a single  point  are  possible.  The  appropriate 
transformation  is  well-known  to  geographers  as  the  Azimuthal 
Equidistant  Projection. 


As  will  be  seen,  a fundamental  step  in  the  algorithm 
is  based  upon  an  iterative  technique  for  the  Euclidean  norm 


k 
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problem  developed  independently  by  Cooper  (1963),  Kuhn  and 
Kuenne  (1962) , and  Weiszfeld  (1936) . A modified  algorithm 
is  used;  namely,  the  Hyperboloid  Approximation  Procedure 
(HAP)  due  to  Eyster  et  al.  (1973),  which  circumvents  the 
difficulties  of  non-differentiability  at  demand  points. 

In  the  plane,  the  original  iterative  procedure  for  the 
single  facility  problem  is  guaranteed  to  converge  to  the 
optimum  location.  A discussion  of  convergence  properties 
is  given  by  Weiszfeld  (1936),  Katz  (1969,  1974)  and  Kuhn 
(1973).  The  iterative  procedure  has  been  found  by  Cooper 
(1963)  and  Eyster  et  al.  (1973),  among  others,  to  solve 
the  problem  with  excellent  results.  With  this  in  mind, 
it  would  seem  advantageous  to  employ  the  technique  in 
solving  the  problem  on  the  sphere. 

The  basic  steps  of  the  approach  follow: 

S P 

(X  on  sphere,  X on  plane) 

STEP  0 : Designate  starting  point  Xq^ . Set  k = 1. 

Set  stopping  criterion  parameter, c. 

STEP  1 : Perform  Azimuthal  Equidistant  Projection 

S 

to  Euclidean  plane  using  X^ as  the  polar  point;  note  that 

c p 

the  point  corresponding  to  X£  ^ is  X^  = (0,0). 

STEP  2 : Employ  HAP  algorithm  using  Euclidean  norm 

P P P 

to  find  global  optimum  Xp,  in  plane.  If  |X^-X^  ^ < e, 

stop.  Otherwise,  go  to  Step  3. 

p 

STEP  3:  Perform  inverse  transformation  on  X,  and 

k 

return  to  sphere  to  get  X^.  If  | <J>(Xj^) -<J>(X^  j)  I < e,  where  <J> 
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is  the  objective  function,  stop.  Otherwise,  let  k = k+1 
return  to  STEP  1. 

In  the  paragraphs  that  follow,  the  mathematical 
bases  for  STEP  1 and  STEP  3 will  be  examined.  STEP  0 
is  discussed  in  Section  III. 12.  For  a detailed  discussion 
of  HAP  (STEP  2)  refer  to  Eyster  et  al.  (1973). 

A specific  advantage  of  this  approach  is  that  while 
solving  the  problem  on  the  plane  one  is  working  with  a 
convex  problem.  Via  the  projection  one  can  consider 
problems  with  demand  points  scattered  around  the  entire 
sphere  without  using  approximations  to  great  circle  dis- 
tances. Also,  under  the  projection  the  number  of  non- 
dif ferentiable  points  is  halved  since  there  is  no  problem 
with  antipodal  points  on  the  plane. 

III.  10.2  Convergence  Properties 

As  stated,  it  is  known  that  the  general  iterative 
technique  converges  to  a global  optimum  on  the  plane. 
Unfortunately,  convergence  of  an  iterative  algorithm  can 
only  be  guaranteed  to  a local  minimum  in  the  general  problem 
on  the  sphere  due  to  its  non-convexity.  Assuming  a con- 
vergent algorithm  is  used  in  the  plane,  this  guarantee  can 
be  made  by  verifying  that  each  iteration  of  the  algorithm 
will  result  in  a strict  improvement  of  the  objective  function 
value.  That  is,  each  time  one  projects  to  the  plane  and 
then  returns  to  the  sphere,  any  movement  in  the  location 
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of  the  servicing  facility  will  result  in  a strict  decrease 
in  total  cost. 

Starting  with  Xq  on  the  sphere  as  a pole,  project 
to  the  plane  using  the  Azimuthal  Equidistant  Projection. 
Suppose  X*  is  the  global  optimum  on  the  plane  for  the 
projected  demand  points.  Assume  it  is  distinct  from  the 
image  of  XQ , X^ . That  is,  assume  Xp  is  not  globally 
optimal  on  the  plane.  If  d^  is  the  image  of  demand  point 
for  each  i and  Xg  is  the  inverse  image  of  the  planar 
optimum  X*  on  the  sphere,  it  follows  from  the  cosine  ine- 
quality for  elliptic  geometry  (Kay  1969)  that  (see  Figure 
3.8)  : 


X*-dil£2  > 


V.,  i = 1 M 


i'  £ 


Sphere 


P lane 


Figure  3.8.  Cosine  inequality  for  elliptic 


geometry . 


All  geodesic  distances  and  angles  from  X are  iso- 

p 


metric  to  the  corresponding  great  circle  distances  and 
angles  on  the  sphere.  All  other  distances  are  distorted  on 
the  plane  to  larger  than  the  actual  distance. 


M M 

Thus  E w.p(Xn,P.)  > E p(X  ,P.)  and  strict 
i=l  1 U 1 i=l  s 1 

improvement  is  guaranteed. 

The  question  arises  as  to  when  one  can  expect  the 
location  of  the  servicing  facility  to  move  during  the  search 
for  a global  optimum  on  the  plane.  For,  if  no  movement  is 
made  from  the  polar  point  X , there  will  be  no  change  in 
the  objective  function  value.  With  this  in  mind,  it  is 
shown  that  the  characteristics  of  the  first  partial  deriva- 
tives on  the  sphere  and  on  the  plane  are  similar.  That  is, 
if  the  partial  derivatives  are  non-zero  in  the  sphere 
problem,  they  will  likewise  be  non-zero  in  the  planar 
problem.  By  employing  a modified  gradient  approach  in  the 
spirit  of  Kuhn  (1973),  there  is  no  difficulty  in  assuming 
differentiability  even  in  the  situation  when  the  polar 
point  is  a demand  point  or  its  antipodal  point. 

Now,  the  convergent  algorithm  in  the  plane  will 
result  in  movement  from  a non-stationary  point  to  a sta- 
tionary point.  So  if  the  foregoing  is  established  one  is 
guaranteed  of  achieving  a strict  improvement  in  any  iteration 
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which  projects  from  a non-stationary  point  on  the  sphere. 

Theorem  3,10.1  Given  an  Azimuthal  Equidistant 

Projection  from  a point  Xq  on  the  sphere,  let  X^  be  the 

corresponding  point  on  the  plane  via  the  projection. 

Suppose  the  first  partial  derivatives  of  the  objective 

function  on  the  sphere  exist  and  are  non-zero  at  Xq . Then 

the  first  partial  derivatives  of  the  planar  objective 

function  are  likewise  non-zero  at  X . 

P 

Proof : The  approach  is  via  contradiction.  Without 

loss  of  generality,  the  coordinates  on  the  sphere  can  be 
translated  so  that  X = (<J>  ,A  ) = (0,0)  in  colatitude, 

S XX 

longitude.  Under  the  projection,  of  course,  X = (X  ,Y  ) = 

p XX 

(0,0)  . 

Suppose  the  first  partial  derivatives  at  X^  vanish 
in  the  plane.  The  proof  will  be  established  by  showing  that 
the  first  partial  derivatives  vanish  for  the  corresponding 
point  Xg  = (0,0)  on  the  sphere. 
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Now,  under  the  projection  it  is  seen  that  for  a 
point  on  the  sphere,  the  corresponding  point 

in  Cartesian  coordinates  is  (a.,b.)  where  a.  = 4>.cosA., 
t>i  = j^sinA^.  This,  of  course,  is  contingent  upon  the 
projection  being  from  = (0,0). 

By  substitution: 


M w.a.  M w . $ . cos  A . m 

£ _ £ i i i _ £ y 

i=l  (a?+b.)2  i=l  [ ((fi.cosA.)^+(<}).sinA.)^]^  i=l  ^ 

li  l l ' i l 


w . <P  . cos  A . 
li  l 


y = Z W. COSA . 


Similarly,  Z w.sinA.  = 0 
i=l  1 

Using  the  form  of  the  objective  function  as  (3.2.7), 

the  following  necessary  conditions  arise: 

a . M -w . (~cos<p  . sinj>+sin<f> . cos<pcos  ( A- A . ) ) 

" i=i  r ================ 

"v  1-  (cost^cosifj+sincf^sin^cos  ( A - A ^ ) ) 

a,(  M -w.  (sin<f»  • sinj)  t sinAcos  A . -cos  AsinA  . ] 

3X  ^ ^ ■ --  ~ 7 

■v  l-(cos<J>^cos <p+s in<jKsin<J>cos(A-Aj.)) 

Looking  specifically  at  when  (<J>,A)  = (<J>  , A ) 

X X 

= (0,0)  it  is  seen  that: 


(0,0) 


M -w^  (cos<Jk  sinO+sin^cosOcos  (-  A^)  ) 
*"  ^ -/rT(cos4»icosO+sin4)isinOcos(-Ai))- 


M -w . sin4>  • cos  A . M 
= Z — -1— ■■  - — — - = Z - w.cosA.  = 0 

1=1  ^ l-cos2<li  1=1 


I 


Also , 
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9 ip 

JX 


(0,0) 


-w.  (sin<l>  .sinOl  sinOcosA  . - cos Osin'*'  . 1) 
E — i - i_ - L_ 

i=l 


= 0 


vn 


^7 

cos  <p , 


This  completes  the  proof. 

The  question  of  whether  the  algorithm  converges 
when  used  with  a convergent  algorithm  in  the  plane  is  not 
really  appropriate  in  this  case  since  the  planar  search  is 
done  with  an  efficient  algorithm  that  uses  a heuristic 
as  a stopping  criterion. 


III. 10. 3 The  Azimuthal  Equidistant  Projection 

In  this  section  the  mathematical  formulations  used 
for  projecting  to  the  planar  surface  are  developed.  This 
is  essentially  STEP  1 of  the  three  major  steps  outlined 
in  Section  III. 10.1. 

Letting  ( <f> , A ) represent  the  spherical  coordinates 
in  colatitude  (0£4>7)  and  longitude  ( - tt < A <_tt ) t the  projection 
from  the  North  Pole  to  the  plane  is  trivial.  Using  polar 
coordinates,  (r,<>),  the  mapping  transformation  is  simply: 

r = <P 
0 = A 

where  0°  longitude,  the  Greenwich  Meridian,  transforms  to 
the  positive  X axis  in  a Cartesian  system,  and  the  South  Pole 
is  a singular  point  arbitrarily  placed  at  7,0).  The 
Cartesian  representation  (x,y)  is  then: 


L 
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x = bcosA 
y = cJ>sinA 

As  mentioned  earlier,  this  particular  projection 
preserves  actual  distances  and  bearing  from  the  North  Pole. 
Error  increases  as  one  determines  distances  between  points 
non-collinear  with  a ray  from  the  origin  and  some  distance 
from  the  pole.  Recognizing  this,  the  algorithm  incorporates 
a bound  which  stops  the  planar  search  whenever  the  movement 
exceeds  a fixed  distance  from  the  pole. 

III. 10.4  Transformation  of  Poles 

In  the  Planar  Projection  Algorithm  (PPA)  one  is 
unfortunately  hardly  ever  projecting  from  the  North  Pole. 

It  is  desired  to  project  from  points  all  over  the  sphere, 
namely,  the  best  point  achieved  in  the  previous  iteration. 
This  requires  a method  to  transform  the  spherical  coordinate 
system  so  that  the  "North  Pole"  is  at  the  desired  pole  for 
projection.  The  transformation  to  the  plane,  then,  is 
trivial . 

Using  a concept  familiar  to  navigators,  it  is  possi- 
ble to  transform  the  coordinate  system  without  ambiguity 
from  one  pole  to  another.  The  concept  to  be  utilized  is 
the  haversine,  where 

hav0  = %(l-cosO)  = sin^(0/2) 

Transferring  from  the  North  Pole  to  new  pole  t*ie 

new  coordinates  (Z^,a^)  for  each  demand  point  P^  = ( 4>  ^ , A ) 
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are  given  by  Haling  (1973) : 

havZ^  = hav|  -i))q | + sin<t>^sin‘t>ghav|  A^-A  | 

havot^  = (havfj -havj  'J)g-Z^  | ) csoJ’qCScZ^ 

Using  inverse  trigonometric  functions,  one  can 
solve  for  (Z^,a^),  the  bearing  and  distance  coordinates 
from  the  new  pole  to  all  demand  points  P^ , i=l,...,M. 

The  end  result,  after  projection  to  the  plane,  is 
a planar  representation  of  the  sphere  with  true  distances 
from  the  new  pole  to  all  points  and  only  "small"  error 
in  distances  from  other  locations  to  the  in  regions  near 
the  pole. 

Solving  for  Z^ : 

havZ^  = hav|  <Jk -4>q  | + sirnfK  sin4>ghav|  A -A  | 

O $ ” f n o A « — A . 

^(l-cosZ/)  = sin  ( — j- ~)  + sin^simlgsin  ( — 

Thus  : 

2 o xn_Ai 

Z^  = Arccos  < 1- 2 (s  in  ( — ~)  + sin^sinij'gsint — ^ — )H  (3.10.1) 
Solving  for  ou  : 

hav«^  = (hav4>^  - hav  | f g-Z^  | ) csccJigCscZ^ 

2 i'i  2 ^o”^i 

^(1-coscu)  = isin  (j- ) - sin  ( — — ) ] csc^gCscZ^ 

Thus  : 

= Arccos{l-2(sinz(^)  - sinz(  - ^ L)  ] cscflgCscZi } (3.10.2) 
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I 

I 
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Although  is  undefined  when  either  <J>q  or  are 
0 or  j this  is  easily  resolved  by  arbitrarily  defining  a. 
to  be  0 since  one  is  at  the  North  or  South  Pole.  is 

defined  as  either  0 or  tt  , whichever  is  correct. 

Observing  that  the  range  of  is  -tt  < ou  £ tt  ( 
it  is  obvious  that  there  will  be  ambiguity  as  to  the  proper 
sign  due  to  the  range  of  [0,it]  for  the  Arccos  function. 

Now,  the  azimuth  a.  is  in  reference  to  the  shortest  great 
circle  arc  from  the  new  pole  to  the  old  pole,  which  is  the 
new  zero  meridian.  In  order  to  determine  the  proper  sign 
it  is  necessary  to  check  the  sign  of  (Aq-A) . The  cases  in 
Table  3. 1 apply. 


| A o - A | <_  ir 

|Aq-A|  > TT 

>- 

0 

1 

>- 

1 V 

o 

— 

a . = a . 
i i 

a . = - a . 
X 1 

Aq-A  < 0 

a.  = -a. 
l l 

a . = a . 

l l 

Table  3.1.  Sign  of  Azimuth. 


An  example  of  the  case  where  |A„-A|  £ tt  and 


Aq-A  ^ 0 is  seen  in  Figure  3.9. 
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Figure  3.9.  Determination  of  Azimuth  (longitude) 

sign. 

111. 10. 5 Returning  to  the  Sphere 

Consider  now  STEP  3 of  the  three  major  steps  in 
Section  III. 10.1.  Prior  to  beginning  the  next  iteration 
and  after  optimization  on  the  plane  it  is  necessary  to 
return  to  the  sphere.  This  is  done  by  converting  the 
globally  optima]  point  (X*,Y*)  to  spherical  coordinates 
us  ing 

<p  = [ (X*)2+(Y*)2j%  (3.10.3) 

A = Arctan(Y*/X*) 

Ambiguity  regarding  the  sign  of  A is  easily  handled.  Once 
on  the  sphere  it  is  simple  to  transform  the  pole  back  to 
the  North  Pole,  using  the  procedure  of  STEP  1.  In  this  way 
the  point  (X*,Y*)  is  expressed  in  spherical  coordinates 
with  the  North  Pole  as  pole . 

1 1 1 . 1 0 . 6 The  Algorithmic  Procedure 

Utilizing  the  results  and  mathematical  development 
of  the  previous  section,  the  following  algorithm  is  obtained. 
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EP1 


EP2 


Terminology  and  notation  used  in  any  iteration,  ITER,  are: 
BOUND  = upper  bound  for  movement  of  search  on  plane 
before  returning  to  sphere 
= stopping  criterion  for  movement  on  the  plane  in 
one  iteration 

= stopping  criterion  for  change  in  objective 
function  value  on  the  sphere  using  (EP2)*F, 
where  F is  the  current  objective  function  value 
= (0,1).  K = 0 if  starting  solution  is  to  be 
projected  centroid.  K = 1 otherwise. 

= bound  on  major  iterations 
= number  of  demand  points 
= weight  of  demand  point  j 
= colatitude  of  demand  point  j in  radians 
= longitude  of  demand  point  j in  radians 


K 

N1 

M 


W. 

J 

♦j 


A . 
J 


1.  Input,  parameters  BOUND,  EP1,  EP2  , K,  N1  , M 
and  demand  points  (b^.A^)  with  weights  W.  for  j = 1,...,M. 

2.  If  K = 0 go  to  step  3.  Otherwise,  go  to  step  4. 

3.  Calculate  the  projected  centroid  (x(|>q,x\q) 
for  starting  solution.  Go  lo  step  5. 

4.  Input  (x^q.xAq). 

5.  Using  (3.2.6)  calculate  initial  objective 
function  value  on  the  sphere  letting  (x1>q,xAq)  be  the 
servicing  facility.  Call  it  OBSPH. 

6.  Set  ITER  = 1,  x<J>  = , and  xA  = xA^. 

7.  If  ITER  _>  Nl,  go  to  step  18.  Otherwise,  go 
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to  step  8. 

8.  Transform  pole  to  (x<!>,xA).  Designate  new 
coordinates  of  demand  points  as  (j’T^.AT^)  for  j = 

Coordinates  of  North  Pole  are  (x<t>^p  > x^p)  ' 

9.  Set  Fb  = ORSPH. 

10.  Transform  all  demand  points  (:T;,'T.)  from 

spherical  coordinates  to  Cartesian  coordinates  via  the 

Azimuthal  Equidistant  Projection.  Designate  them  as 

(A. , B . ) . 

J J 

11.  Using  the  Hyperboloid  Approximation  Procedure 

(HAP) , or  any  other  Euclidean  norm  single  facility  algorithm, 
solve  the  problem  in  E“  for  the  global  minimum  (X*,Y*) . 

12.  Check  to  see  if  the  first  stopping  criterion 

is  satisfied.  If  [(X*)^+(Y  )^]2  < EP1  then  do  steps  13 

through  15  and  go  to  step  18.  Otherwise,  go  to  step  13. 

13.  Transform  (X*,Y*)  from  Cartesian  to  spherical 
coordinates,  thus  returning  to  the  sphere  with  (x<f>,xA) 

as  pole  to  get  (x<j'*,xA*). 

14.  Transform  the  pole  at  (x<J>,xA)  back  to  the  North 
Pole  (xj^p , xA^p)  to  get  (x<|>* , xA*)  in  terms  of  the  North 
Pole.  Call  the  new  point  (x<Ji,xA). 

15.  Using  (3.2.6)  evaluate  the  objective  function 
OBSPH  on  the  sphere,  letting  (xj’.xA)  be  the  servicing 
facility.  Go  to  step  16  (unless  stopping  criterion  in 
step  12  is  satisfied) . 


16, 


5et  Fa  = 


OBSPH. 
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17.  Check  to  see  if  the  second  stopping  criterion 
is  satisfied.  If  | FA - Ffi | < (EP2)*FA>  go  to  step  18. 
Otherwise,  go  to  step  7. 

18.  Stop. 

III. 11  A Cyclic  Meridian  Search  Algorithm  (CMS) 

111. 11.1  Introduction 

Recognizing  that  an  appreciable  amount  of  time  is 
required  for  calculation  of  the  gradient  in  the  planar 
search  and  in  transforming  back  and  forth  from  the  sphere, 
one  can  reasonably  consider  the  possibility  of  both  search- 
ing entirely  on  the  spherical  surface  and  avoiding  the 
use  of  derivatives.  Another  factor  to  encourage  such  an 
approach  is  that  due  to  the  way  a gradient  search  functions , 
it  will  proceed  downhill  to  a local  optimum  once  it  is 
within  its  "region  of  attractiveness."  A derivative-free 
approach  may  be  able  to  avoid  local  optima. 

Two  questions  come  immediately  to  mind.  That  is, 
how  should  one  choose  a search  direction,  and  how  should 
the  line  search  be  made.  Two  simple  procedures  were  used. 

1 1 1 . 1 1 . 2 The  Search  Direction 

Concerning  search  directions,  the  developed 
algorithm  always  searches,  in  each  cycle,  in  a direction 
orthogonal  to  a meridian  along  a great  circle  track  and 
then  along  a meridian.  This  cycle  is  repeated  until  an 
initial  stopping  criterion  is  satisfied.  At  that  juncture 


62 


random  search  directions  are  tried  in  order  to  further 
improve  the  objective  function  value.  This  procedure 
helps  to  avoid  getting  stranded  on  ridges  and  in  valleys- - 
as  can  happen  in  a cyclic  coordinate  (one  at  a time) 
search  on  the  plane. 

The  occurrence  of  such  ridges  and  valleys  on  the 
sphere  can  best  be  demonstrated  via  example.  Consider 
the  four  point  single  facility  location  problem  using  data 
set  D4  of  Appendix  B,  depicted  in  Figure  3.10.  The  global 
optimum  in  degrees  of  latitude  and  longitude  is  at  approxi- 
mately (33,57),  yet  a cyclic  search  that  does  not  incorporate 
random  search  directions  could  get  stuck  at  the  point  (0,20) 
even  though  it  is  clearly  not  a local  minimum.  This  would 
happen  if  the  optimization  started  with  20°  longitude  as  a 
fixed  meridian  and  a line  search  along  it  converged  to 
(0,20).  At  this  point  an  orthogonal  search  would  be  along 
the  great  circle  arc  which  is  identical  to  the  Equator, 
or  0°  latitude.  Such  a search  would  not  decrease  the  value 
of  the  objective  function  so  the  algorithm  would  terminate 
at  (0,20). 

The  cyclic  search  process  can  be  quite  slow  around 
ridges  and  valleys  due  to  a tendency  of  the  search  to  zig-zag 
along  them.  This  is  discussed  with  regard  to  this  four  point 
problem  in  Chapter  V.  Overall,  however,  the  algorithm 
seemed  to  be  quite  efficient  when  compared  to  the  planar 
projection  method. 


It  is  easy  to  visualize  that  although  search  direc- 
tions are  orthogonal  within  each  cycle,  the  actual  direction 
of  search  on  the  sphere  can  be  in  any  direction- -except  in 
a pathological  case  when  the  line  search  in  a fixed  direc- 
tion moves  an  exact  multiple  of  tt / 2 at  every  iteration. 

This  property  in  itself  would  appear  to  make  the  technique 
more  powerful  than  the  counterpart  cyclic  coordinate 
descent  algorithm  on  the  plane  (Zangwill  1974) . 

1 1 1 . 1 1 . 3 The  Line  Search 

Regarding  the  line  search,  it  must  be  recognized 
that  although  the  search  along  a geodesic  (meridian)  in 
any  direction  from  a pole  is  bounded  by  it,  the  objective 
function  is  not  in  general  unimodal  over  the  search  region. 
As  a result,  the  standard  efficient  line  search  techniques 
(Fibonnaci,  Golden  Section)  will  not  guarantee  a global 
optimum.  Also,  since  local  optima  exist  in  the  general 
problem,  one  should  give  thought  to  techniques  of  escaping 
them  when  at  all  possible.  It  turns  out  that  an  adaptation 
of  a simple  search  technique  due  to  Bazaraa  (1975)  will 
perform  a good  line  search  and  in  many  cases  avoid  local 
optima,  especially  when  the  local  optima's  region  of  attrac- 
tion are  relatively  small  and  the  objective  function's 
surface  is  not  too  flaL. 

The  algorithmic  steps  of  the  line  search  will  be 
described  in  the  next  section.  A general  outline  follows: 
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Given  the  objective  function  where  f(<J>*,A)  = Z 
where  A is  a fixed  meridian--along  which  the  line  search 
is  being  made--check  to  see  if  f($*+e,A)  or  f(<J>*-e,A) 
results  in  a decrease  in  the  objective  function  value  for 
arbitrarily  small  c > 0.  If  neither  decrease  the  objective 
function  value  tne  search  is  terminated.  Suppose  it  is 
determined  that  f(6*+e,A)  is  the  appropriate  search  direc- 
tion. The  intent  of  the  search  then  is  to  increase  the 
colatitude  $ as  much  as  possible  for  fixed  A,  all  the  while 
decreasing  the  objective  function  value.  Choosing  an 
initial  step  size  A > 0,  and  an  initial  acceleration 
factor  S = 1,  evaluate  f (4>*+SA , A)  . If  f(<j>*+SA,A)  < Z* , 
then  Z*  is  replaced  by  f(4>*+SA,A),  <}>*  is  replaced  by  <|)*+SA, 
S is  replaced  by  uS  (a>l) , and  the  process  is  repeated 
until  the  first  failure  is  encountered.  If  only  f ($*-£, A) 
had  decreased  the  objective  function  value,  then  S would  be 
replaced  by  -1  and  the  colatitude  41  would  be  decrea.  ed  as 
much  as  possible. 

When  the  first  failure  is  encountered,  the  step 
size  A is  decreased  and  the  process  continued.  This 
reduction  in  step  size  takes  place  after  each  failure  until 
the  minimum  step  size  results  in  a failure.  In  this  case, 
an  approximate  minimum  (possible  local)  of  f is  4>* . 

A prominent  feature  of  this  approach  is  the  fact 
that  the  initial  step  size  A is  a parameter  that  can  be 
controlled  so  that  in  many  cases  it  is  possible  to  step 
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out  of  the  range  of  a local  optimum.  A graphical  picture 
of  this  is  in  Figure  3.11.  Note  that  the  first  evaluation 
gives  f(<j>*+A)  < f (4>*)  and  the  local  minimum  is  escaped. 

It  must  be  emphasized  that  selection  of  A has  a pronounced 
effect  on  success  in  escaping  local  minima,  and  little  is 
known  a priori  as  to  appropriate  choices. 


Figure  3.11.  Escaping  a local  minimum. 

I I I . 1 1 . 4 The  Algorithmic  Procedure 

Armed  with  the  foregoing  preliminary  comments,  the 
Cyclic  Meridian  Search  (CMS)  Algorithm  is  now  described. 
The  following  notation  and  terminology  is  used: 

S = initial  acceleration  factor  (=1) 

a = acceleration  factor  for  step  size  (>1) 

N = switching  factor  for  reduction  of  basic  step 


size 
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e 


A 1 

AMIN 
A INC 
K 

M 

W. 

J 

♦j 

A . 

J 

ITBOUN 

EPI 

EP2 

LRANBO 


= small  seep  size  used  to  find  direction  (±)  of 
search  on  fixed  meridian 
= initial  step  size 
= minimum  step  size 
= incremental  step  size 

= (0,1),  K=0  if  projected  centroid  is  to  be 
starting  solution;  K=1  otherwise. 

= number  of  demand  facilities 
= non-negative  weight  associated  with  demand 
point  j 

= colatitude  of  demand  point  j in  radians 
= longitude  of  demand  point  j in  radians 
= bound  on  number  of  cyclic  searches 
= stopping  criterion  for  search  movement 
= stopping  criterion  for  improvement  in  objective 
function  using  EP2*(F)  where  F is  the  current 
objective  function  value 
= number  of  random  search  directions  employed 
once  either  stopping  criterion  is  satisfied 


1.  Input  parameters  a,  M,  Al,  AMIN,  AINC,  K,  c, 

ITBOUN,  EPI,  EP2,  LRANBO,  and  demand  points  ($.,A  ) with 
weights  W.  for  j = 1 M. 

2.  If  K = 0,  go  to  step  3.  Otherwise  go  to  step  4. 

3.  Calculate  the  projected  centroid  (x<1q,xAq)  for 
starting  point.  Go  to  step  5. 

4.  Input  starting  point  (x^q.xAq). 
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5.  Set  <t>*  = x<l>0>  = xAQ' 

6.  Calculate  initial  objective  function  value  on 
the  sphere  using  (3.2.6),  letting  (<J>*,A*)  be  the  servicing 
facility.  Call  it  F^. 

7.  Fix  the  meridian  at  A*,  and  perform  a line 
search  from  <1*  to  get  (x4>,xA)  and  F , the  objective  function 
value  (an  improvement  over  F j ^ ) . 


8. 

Set 

FIN  = Fs‘ 

9 . 

Set 

ITER  = 0. 

10. 

Set 

ITER  = ITER 

+ 1 , LRAN  = 0 . 

11. 

Set 

up  a search 

from  the  current  pole  along 

the  meridian  i n / 2 . Set  (ij>*,A*)  = ( 0 , tt  / 2 ) . 

12.  See  if  the  bound  on  iterations  has  been  exceeded. 
If  ITER  is  less  than  ITBOUN,  go  to  step  13.  Otherwise  go 

to  step  33. 

13.  Transform  the  pole  to  (x<J>,xA)  using  procedures 
in  Section  III..10.4  Store  North  Pole  as 

(1*’np,anp^  is  re<luired  in  stGP  28- 

14.  Fixing  the  meridian  at  A*,  perform  a line 

search  from  <4*  to  get  (x<t>,xA)  and  , the  corresponding 
objective  function  value. 

15.  Check  to  see  whether  a random  search  direction 
has  been  used  within  the  iteration.  If  LRAN  is  greater 
than  zero , go  to  step  18.  Otherwise  go  to  step  16. 

16.  Insure  that  two  orthogonal  searches  are  per- 
formed before  checking  stopping  criteria.  The  two  searches 
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complete  a cycle  in  the  two  dimensional  spherical  space. 

In  each  cycle  the  first  search  is  performed  along  a great 

circle  arc  orthogonal  to  a meridian  when  the  North  Pole 

is  the  polar  point,  and  the  second  search  is  along  a 

meridian.  The  pattern  of  such  a search  is  depicted  in 

Figure  3.12.  Observe  that  as  the  region  of  search  gets 

smaller,  the  search  closely  approximates  a cyclic  coordi- 
2 

nate  search  in  E . 

If  X*  - n/2  go  to  step  17.  Otherwise  go  to  step  21. 


Figure  3.12.  Sample  search  pattern  for  cyclic 
meridian  search. 


17.  Set  up  values  for  use  in  testing  stopping 
criteria  after  completion  of  second  half  ol  cycle  Set 
Fj  = F^  - F^j , x<(>  ^ = x-i>  and  - F^ . Go  to  step  28. 

18.  If  a random  search  direction  has  been  used 
within  the  iteration,  a stopping  criterion  is  tested  to 
see  if  it  is  still  satisfied.  If  jF^-F^I  _<  EP.  (F^)  go  to 
step  22.  Otherwise  go  to  step  19. 
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20.  The  second  stopping  criterion  is  tested  after 
use  of  a random  search  direction  to  see  if  it  is  still 
satisfied.  If  |xif>|  £ EP1  go  to  step  25.  Otherwise  go 

to  step  28,  and  begin  a new  iteration. 

21.  This  step  is  reached  after  one  cycle  of  the 
line  search  is  completed- -namely  two  orthogonal  searches. 

A check  is  made  to  see  if  a stopping  criterion  is  satisfied. 
If  I F^“+(F  -F  ) " ' 2 EP2(F^)  go  to  step  22.  Otherwise  go 

to  step  23. 

22.  Set  Fg  = F^  and  go  to  step  25. 

23.  Set  F = FM. 

s N 

24.  The  second  stopping  criterion  is  checked  after 

2 2 k 

a cycle  of  the  line  search  is  completed.  If  (xf^+xf'  ) 2 _< 
EP1  go  to  step  25.  Otherwise  go  to  step  28  and  begin  a 
new  iteration. 

25.  This  step  is  reached  whenever  any  of  the  stop- 
ping criteria  are  satisfied.  It  enables  the  algorithm  to 
make  a search  in  a random  direction.  Set  LRAN  = LRAN  + 1. 

26.  Check  to  see  if  the  upper  bound  on  random 
search  directions  in  any  iteration  has  been  exceeded. 

If  LRAN  is  greater  than  LRANBO  go  to  step  28.  Otherwise 
go  to  step  27. 

27.  Generate  a random  number  RAN,  -1  £ RAN  £ 1, 
set  t>*  = 0 and  let  the  fixed  meridian  (line  of  search) 
be  A*  = RAN(h).  Go  to  step  13. 


28.  Convert  (xif>,x\)  to  coordinates  with  the 
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North  Pole  as  pole  using  ($^p,^„p)  from  step  13  and  pro- 
cedures of  Section  III.10A  for  transforming  poles.  This  is 
done  in  preparation  for  performing  the  next  iteration  or 
cycle . 

29.  Evaluate  the  objective  function  value  on  the 
sphere  using  (3.2.6)  and  letting  (x;,xO  be  the  servicing 
facility . 

30.  If  LRAN  is  greater  than  LRANBO,  go  to  step  33. 
Otherwise  go  to  step  31. 

31.  Determine  whether  the  line  search  cycle  has 
been  completed.  If  \*  = n/2,  go  to  step  32.  Otherwise  go 
to  step  10. 

32.  Begin  the  second  half  of  the  line  search 
cycle.  Set  (<fi*,A*)  = (0,0)  and  go  to  step  13. 

33.  Stop. 

The  line  search  procedure,  being  the  heart  of  the 
algorithm,  is  now  described  as  a separate  entity. 

1.  Input  parameters  M,  A1  , AMIN,  AINC,  cx , N, 
(<t>*,A*),  > and  demand  points  (<J>.,Aj)  with  weights  W.  . 

Let  Fin  = f(<t>*,A*). 

2.  Set  ICOUNT  = 0.  This  counter  is  used  in 
determining  whether  the  appropriate  direction  of  search 
on  the  fixed  meridian  has  been  determined. 

3.  Set  R = 1.  R is  a parameter  which  fixes  S 

at  a proper  value,  t 1,  after  each  change  in  step  size  A. 

4.  Set  A = A1 . 


72 


5.  Set  S = R. 

6.  Set  Kq  = 0 . Kq  counts  the  number  of  accelera- 
tions for  a given  step  size. 

7.  K,  = Kq  + 1. 

8.  Check  to  see  how  many  steps  have  been  made  for 
a particular  step  size  A.  This  permits  one  to  restrict 
the  search  to  a maximum  size  interval.  For  example,  if 

Al  = .1  radians,  S = 1,  a = 2 and  N = 5,  the  largest  move 
along  a line  due  to  acceleration  would  be  2^(.l)  = 3.2, 
which  is  just  greater  than  v.  If  Kq  is  greater  than  N, 
go  to  step  36.  Otherwise  go  to  step  9. 

9.  If  ICOUNT  = 0,  go  to  step  27.  Otherwise  go 
to  step  10. 

10.  Set  ICOUNT  = ICOUNT  + 1. 

11.  Set  X'i^  = <J>*+S(A)  and  xA  = A*. 

12.  Due  to  the  bound  on  colatitude,  it  is  neces- 
sary to  insure  that  x4>^  is  never  less  than  zero  or  greater 
than  it.  If  xif)„  is  less  than  or  equal  to  it  and  greater  than 
or  equal  to  zero,  go  to  step  28.  Otherwise,  go  to  step  13. 

Steps  13  through  22  re-orient  the  direction  of  search  within 
the  bounds  of  colatitude  and  along  the  proper  meridian  when 
necessary . 

13.  If  x>(>^  is  less  than  or  equal  to  0 , go  to  step 
Otherwise  go  to  step  14. 

14  Set  xfN  = 2tj-x4>n. 
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15.  If  A*  + v is  greater  than  it,  go  to  step  17. 


Otherwise  go  to  step  16. 


16. 

17. 

18. 

19. 

Otherwise  go 

20. 


Set  xA  = A*  + v and  go  to  step  22. 

Set  xA  = A*  - it  and  go  to  step  22. 

Set  x*Dn  = 

If  A*  f 7i  is  greater  than  , go  to  step  21. 
to  step  20. 

Set  xA  = AT?  + 77  and  go  to  step  22. 


21.  Set  xA  = A*  _ A. 

22.  Set  R = -R  and  S = -S . 

23.  Evaluate  the  objective  function  at  (x4>^,xA) 
using  (3.2.6).  Call  it  F. 

24.  If  ICOUNT  =0,  go  to  step  31.  This  means  that 
the  appropriate  direction  for  search  along  xA  has  not  yet 
been  determined.  Otherwise,  go  to  step  25. 

25.  Check  to  see  if  a decrease  in  the  objective 
function  value  has  been  obtained.  If  not,  then  the  step 
size  must  be  reduced.  If  (F  - F)  is  less  than  zero  go  to 
step  26.  Otherwise  go  to  step  33. 

26.  Set  R = -R  and  S -S , and  go  to  step  37. 

27  Find  the  appropriate  direction  (t)  of  search 
from  ji*  by  testing  the  effect  on  the  objective  function  due 
to  a small  movement  from  <i>*  to  4>*  + e along  A*.  Set  x<  ^ 

= <J>*  + t and  xA  = A*.  Go  to  step  12. 

28.  Evaluate  the  objective  function  at  (x^^.xA) 
using  (3.2.6).  Call  it  F. 
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29. 

to  step  30. 

30. 

decrease  in 

(F  - F)  is 
s 

to  step  33. 


If  ICOUNT  = 0 , go  to  step  31.  Otherwise  go 

Check  the  stopping  criterion  to  see  whether  a 
the  objective  function  has  been  obtained.  If 
less  than  zero,  go  to  step  36.  Otherwise  go 


the 

zero 


F = 
s 


31.  This  step  is  entered  only  when  determining 
initial  direction  of  search.  If  (F  - F^)  is  less  than 
, go  to  step  33.  Otherwise  go  to  step  32. 

32.  Set  R = -R,  S = -S,  ICOUNT  = ICOUNT  + 1, 

F,^  and  go  to  step  11. 

33.  Set  ICOUNT  = ICOUNT  + 1,  F = F,  41*  = x<|>N , 


A*  = xA . 


34.  Check  to  see  if  the  proper  direction  (f)  has 
just  been  found.  If  ICOUNT  = 1,  go  to  step  10.  Otherwise 
go  to  step  35. 

35.  Set  S = a(S)  and  go  to  step  7. 

3b.  Set  A*  = xA . 

37.  Decrease  step  size  A by  AINC. 

38.  If  A is  less  than  AMIN,  go  to  step  39.  Other- 
wise go  to  step  5. 

39.  Set  X'f>  = <f  ' xA  = A*. 

40 . Stop . 


III. 12  Starting  Solutions 

A fundamental  consideration  in  any  iterative  tech- 
nique is  what  to  use  as  a starting  solution.  In  this 
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section  the  characteristics  of  a number  of  approaches  will 
be  discussed. 

1 1 1 . 1 2 . 1 Projected  Centroid 

In  testing  the  two  algorithms,  it  was  found  that 

the  projected  centroid  provided  a "good"  starting  solution. 

Often  it  was  in  the  vicinity  of  the  global  optimum,  and  in 

many  cases  the  global  optimum  was  achieved  when  this  point 

started  the  iterative  procedure.  Finding  the  centroid 
3 

m E for  the  Euclidean  norm  and  projecting  it  to  the 
sphere's  surface  along  a vector  emanating  from  the  sphere's 
center  is  easy  to  do  analytically.  As  discussed  in 
Section  III.  8,  this  point  was  developed  as  an  approximate 
solution  by  Wendell.  Intuitively  it  would  seem  to  be  a good 
starting  solution  for  a problem  with  demand  points  all  over 
the  sphere  since  the  projected  centroid  would  at  least 
locate  one  in  a hemisphere  containing  a majority  of  the 
weighted  demand. 

3 

1 1 1 . 1 2 . 2 Projected  Optimum  for  Euclidean  Norm  in  E 
Although  comparatively  more  expensive  to  determine, 

this  point  could  conceivably  be  better  than  the  projected 
centroid  in  certain  cases.  However,  it  unfortunately 
requires  an  iterative  technique  to  find  it. 

111. 12. 3 Random  Start 

One  could  generate  random  starting  points  on  the 


spherical  surface.  This  is  the  basis  for  the  most  practical 
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way  to  help  ascertain  whether  global  or  local  optima  have 

been  achieved,  other  than  by  plotting  contours.  In  the 

special  case  when  either  the  centroid  or  Euclidean  distance 
3 

solution  in  E is  at  the  center  of  the  sphere,  this  procedure 
would  be  an  alternative. 

I I I . 1 2 . 4 Demand  Points  and  their  Antipodal  Points 

Starting  the  iterative  procedure  at  a demand  point 
or  its  antipodal  point  can  be  informative.  Since  at  most 
only  a local  optimum  can  be  expected  in  an  iterative  tech- 
nique, a number  of  different  starting  points  will  likely 
be  tried.  Use  of  a demand  point  will  reveal  quickly  whe- 
ther it  is  a local  optimum  or  not.  If  it  is,  future 
starting  points  should  be  selected  at  some  distance  from 
it  in  order  to  try  to  avoid  converging  to  it  again.  The 
antipodal  point  certainly  fits  this  criterion. 

I I I . 1 3 Avoiding  Local  Minima 

It  is  well  known  that  a guarantee  of  an  obtainable 
global  optimum  occurs  only  when  the  search  region  is  convex 
and  when  the  objective  function  is  unimodal  in  the  appropri- 
ate form,  i.e.,  convex  for  a minimization  problem.  In 
working  with  a non-convex  programming  problem,  as  in  this 
case,  one  must  be  concerned  with  methods  of  avoiding  local 
optima.  Reklaitis  and  Phillips  (1975)  state  that  all  known 
techniques  except  those  employing  statistical  sampling 
techniques  will  generate  only  local  minima.  They  reference 
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work  by  Clough  (1969)  and  Liau  (1973)  concerning  the  sta- 
tistical approach. 

A standard  approach  when  using  iterative  (numerical) 
techniques  is  to  accept  the  fact  that  only  local  minima 
can  be  obtained  and  solve  the  problem  a number  of  times 
using  random  starting  points.  One  could  then  state  with 
some  degree  of  confidence  that  the  global  optimum  has 
been  achieved- - the  best  solution  obtained  thus  far.  How- 
ever, there  is  usually  no  way  of  knowing  how  many  local 
minima  exist  and  the  possibility  of  a global  optimum  lying 
elsewhere  cannot  be  disregarded. 

Intuitively  it  is  seen  that  many  starting  points 
will  result  in  the  same  local  optimum  if  they  are  within 
the  "hollow''  of  that  local  minimum.  Recognition  of  this 
leads  to  another  class  of  approaches.  One  could  use  a 
solution  for  one  problem  and  employ  various  techniques  for 
"jumping  out"  of  the  range  of  the  current  local  minimum. 
Hesse  (1973)  developed  a heuristic  procedure  which  fits  in 
this  class.  A penalty  function  based  upon  an  added  spheri- 
cal constraint  is  used  in  order  to  get  away  from  a local 
minimum.  In  this  way  an  attempt  is  made  to  find  better 
local  optima. 

A third  class  of  methods  attempts  to  gain  informa- 
tion about  the  entire  search  region  (in  the  case  at  hand, 
either  the  spherically  convex  hull  or  the  entire  sphere) . 
Attention  would  gradually  be  concentrated  upon  smaller 
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regions  which  appear,  under  some  criterion,  to  be  likely 
to  contain  the  global  optimum.  An  approach  using  networks 
that  belongs  to  this  class  has  been  suggested  by  Robinson 
(1972).  Hartman  (1973)  examined  six  variations  of  the 
three  classes,  restricting  consideration  to  the  "essentially 
unconstrained"  nonlinear  programming  problem: 

min  f (x) 

subject  to  x c S C En 

where  the  boundaries  of  S do  not  determine  the  solution. 

The  restriction  to  En  is  not  necessary,  so  his  findings 
can  be  applied  to  the  problem  at  hand. 

In  general,  Hartman  found  that  for  the  variations 
used,  the  first  two  classes  performed  better  than  the  third. 
On  difficult  problems,  though,  even  the  best  of  the  methods 
will  frequently  fail  to  locate  the  global  optimum.  Hartman 
also  found  that  the  best  results  were  obtained  by  methods 
which  do  the  least  random  searching. 

The  subject  of  avoiding  local  minima  is  a fertile 
field  for  further  research.  Approaches  for  avoiding  local 
minima  can  be  adapted  successfully  to  capitalize  on  fea- 
tures of  a particular  problem,  and  can  even  be  incorporated 
within  the  search  process  itself.  One  feature  of  the  cyclic 
search  algorithm  of  Section  III. 11  is  its  ability  to  jump 
out  of  the  range  of  a local  minimum  during  a line 
search.  The  success  of  this  "escape  mechanism"  is  a function 
of  the  initial  step  size  (parameter  Al)  and  the  relative 
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deepness  of  the  hollow  containing  the  global  optimum  as 
compared  to  those  containing  other  local  minima. 

III. 14  Application  to  Multifacility  Location-Allocation 
Problems 


III. 14.1  Model  of  Problem 

In  this  section  a specific  class  of  location  prob- 
lem, the  location-allocation  problem,  will  be  considered. 
Surveys  of  work  on  this  problem  have  been  presented  by 
Scott  (1970)  and  Cooper  (1964) , and  an  extensive  bibliography 
has  been  accumulated  by  Lea  (1973) . A general  objective  in 
this  class  of  problem  involves  determination  of  the  number 
of  new  servicing  facilities,  their  location,  and  optimum 
allocation  of  demand  points  to  servicing  facilities.  In 
this  section,  however,  it  will  be  assumed  that  the  number 
of  new  facilities  to  be  located  is  known. 

A mathematical  formulation  of  this  location-allocation 
problem  for  the  Euclidean  norm  is: 
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= weights  of  demand  facilities 
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(a^.bj)  = P.  are  demand  point  locations 
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= number  of  demand  points 
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number  of  facilities  to  be  located 
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As  exhibited  by  Cooper  (1967)  this  problem  in  the 
plane  is  neither  convex  nor  concave,  and  is  hence  much  more 
difficult  to  solve  than  the  convex  single  facility  problem. 
The  multifacility  location-allocation  problem  on  the  sphere 
is  no  less  difficult. 

Two  approaches  to  the  problem  in  the  plane  are 


evident,  and  have  been  developed  over  the  years.  One  is 
via  combinatorial  programming,  recognizing  the  0-1  nature 
of  the  variables  ' ^ . A modification  of  this  approach  to 
solve  the  counterpart  multifacility  problem  on  the  sphere 
is  the  intent  of  this  section.  The  other  major  approach, 
used  by  Cooper  (1967),  is  via  an  extension  of  the  single 
source  algorithm.  Although  such  an  approach  can  only 
obtain  local  minima,  it  is  very  inexpensive  with  regard 
to  computer  time,  and  is  the  only  efficient  approach  when 
solving  very  large  problems. 

Mathematical  formulation  of  the  spherical  multi- 
facility problem  is,  of  course,  similar  to  the  planar 
formulation,  with  the  exception  of  the  different  distance 
measure  and  the  "spherical"  equality  constraint  which 
restricts  solutions  to  the  spherical  surface. 
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where  (a.,  b.,  c.)  = P.  are  demand  points,  and  all  other 
J J J J 

terminology  is  as  in  the  planar  formulation  (3.14.1). 

I I I . 1 4 . 2 A Planar  Algorithm 

Kuenne  and  Soland  (1971,  1972)  provide  a technique 
for  solving  the  location-allocation  problem  via  branch  and 
bound.  A fundamental  aspect  of  the  algorithm  is  sequential 
solution  of  a number  of  single  source  problems,  for  which  it 
is  possible  to  find  a global  optimum  in  the  plane. 

One  can  consider  direct  application  of  a great  circle 
metric  single  facility  algorithm  in  the  branch  and  bound 
algorithm.  However,  since  only  local  optimality  is  obtain- 
able in  the  general  case  for  the  single  facility  problem 
on  the  sphere,  one  can  expect  no  better  than  a local  minimum 
for  the  multifacility  problem.  Also,  as  will  be  seen,  the 
lower  bound  calculated  in  the  branch  and  bound  algorithm 
may  not  be  an  actual  lower  bound.  Modifications  can  be 
made,  though,  to  partially  resolve  these  difficulties. 

It  is  interesting  to  note  that  Kuenne  and  Soland 
recognized  the  possible  impact  of  great  circle  distances 
and  devised  an  approximation  technique  for  problems 
covering  large  regions  of  the  earth.  It  is  based  upon  a 
Mercator  projection,  where  distances  are  computed  as  rhumb- 
line map  distances  and  then  converted  via  an  approximation 
to  great  circle  distance  (by  multiplying  the  distances  by 
the  cosines  of  the  midpoints  of  the  spherical  coordinate 
latitudes) . Kuenne  and  Soland  hold  that  the  approximation 
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will  work  well  for  most  purposes.  The  large  problem 
solved  by  Kuenne  and  Soland  covered  the  contiguous  United 
States,  for  which  a standard  Mercator  projection  is  suitable. 
However,  such  a projection  will  fail  dramatically  if  demand 
points  are  near  the  polar  regions  or  in  regions  larger 
than  a hemisphere. 

III. 14. 3 Modifications  for  Spherical  Problem 

Inherent  in  Kuenne  and  Soland' s algorithm  is  a 
presumption  that  global  optimality  can  be  guaranteed  for 
the  single  facility  problem.  Single  facility  sub-problems 
are  solved  and  the  results  are  used  in  development  of 
lower  bounds,  feasible  solutions  and  branching  points  at 
each  iteration.  Noting  that  the  contiguous  United  States 
can  be  located  within  a disk  of  diameter  rrR/2,  where 
R is  the  earth's  radius,  one  can  be  confident  of  obtaining 
a global  optimum  (upon  accepting  the  validity  of  Conjecture 
3.7.1).  It  is  evident  that  Kuenne  and  Soland  implicitly 
assumed  the  validity  of  the  conjecture. 

Unfortunately,  in  the  general  case  global  optimality 
cannot  be  guaranteed  for  the  single  facility  problem.  This 
leads  to  difficulties  in  determining  a lower  bound  for  the 
branch  and  bound  algorithm.  As  used  by  Kuenne  and  Soland, 
these  bounds  consist  of  the  sum  of: 

1.  cost  contributions  arising  from  optimally 
locating  r servicing  facilities  to  which  more  than  two 
demand  points  have  been  assigned,  and 


2.  underestimates  of  costs  associated  with 


optimally  locating  all  servicing  facilities  to  service  all 
unassigned  demand  points  in  addition  to  all  demand  points 
which  had  facilities  servicing  only  them. 

The  underestimates  of  (2)  present  no  problem  since 
they  are  based  upon  inter-demand  point  distances.  No 
optimization  problem  is  involved.  However,  since  only  a 
local  minimum  can  be  guaranteed  in  general  for  the  location 
of  each  of  the  r servicing  facilities  in  (1) , there  is  no 
way  to  discern  whether  the  sum  of  differential  costs  between 
the  r global  and  local  solutions  exceeds  the  differential 
amount  between  the  actual  costs  of  (2)  and  the  calculated 
underestimate  of  (2) . If  it  does  exceed  the  differential 
amount,  the  calculated  lower  bound  would  not  be  an  actual 
lower  bound.  In  other  words,  the  costs  associated  with  the 
local  solutions  may  be  so  "bad"  as  to  overcome  the  under- 
estimated costs  of  (2) . 

Considering  the  general  case,  then,  it  is  necessary 
to  modify  the  Kuenne  and  Soland  algorithm  in  order  to  make 
it  adaptable  to  the  sphere  problem.  It  is  important  to 
point  out,  however,  that  validity  of  Conjecture  3.7.1 
will  permit  direct  application  of  the  branch  and  bound 
algorithm  (with  the  spherical  single  facility  location 
algorithm  as  a subroutine)  when  all  demand  points  can  be 
contained  within  a disk  of  diameter  £ ttR/2  of  their 
servicing  facility. 
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The  procedural  steps  for  the  general  case  required 
at  each  iteration  follow: 

For  each  tentatively  located  servicing  facility  K, 
check  to  see  if  all  points  allocated  to  K are  contained 
within  a disk  of  diameter  <_  nR/2  (this  determination  of 
the  covering  circle  size  is  accomplished  via  the  pro- 
cedures in  Chapter  IV).  If  so,  solve  for  the  global 
optimum  and  use  the  objective  function  value  for  the 
lower  bound  contribution.  If  not,  then: 

a.  solve  for  the  local  optimum  and  use  the  solution 
point  as  an  input  in  determining  the  branching  point  and 
feasible  solution. 

3 

b.  solve  the  unconstrained  problem  in  E for  the 
Euclidean  norm.  Use  the  resulting  objective  function 
cost  for  a lower  bound  contribution  (see  Section  III. 9). 
As  is  evident,  the  lower  bounds  will  not  be  as  tight  as 
in  the  basic  algorithm,  so  fathoming  will  not  generally 
occur  as  early. 

Note  that  once  a node  is  reached  where  it  cannot 
be  concluded  that  a global  optimum  is  found,  no  further 
partitioning  from  that  node  will  ever  result  in  a global 
optimum.  This  is  because  each  branch  assigns  more  demand 
points  to  a servicing  facility,  not  less.  Hence,  the  disk, 
or  covering  circle,  can  only  increase  in  size. 

The  advantage  of  the  preceding  approach  is  that  a 
very  good  local  minimum  will  eventually  result  because  of 
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the  implicit  examination  of  all  possible  demand  point- 
servicing facility  allocations. 

Even  though  incorporation  of  the  single  facility 
algorithm  for  the  sphere  as  a basis  for  the  Kuenne  and 
Soland  branch  and  bound  approach  requires  weakening  of  the 
lower  bounds,  it  is  much  better  than  the  Mercator  approxi- 
mation for  a problem  covering  a region  as  large  as  the 
Eurasian-African  land  mass.  As  already  mentioned,  distances 
between  extreme  demand  points  would  be  distorted  significantly 
in  a Mercator  Projection  using  rhumb  lines  as  an  approxi- 
mation to  great  circle  distances.  The  single  source 
algorithms  which  have  been  presented  are  not  affected  by 
the  size  of  the  region  in  determination  of  inter- demand 
point  distances,  nor  are  they  affected  by  the  region's 
location  on  the  sphere. 

I I I . 1 A . 4 Another  Approach 

As  with  any  branch  and  bound  algorithm,  large 
scale  problems  are  generally  not  able  to  be  efficiently 
solved.  An  increase  in  servicing  facilities  (N)  along 
with  a corresponding  increase  in  demand  points  (M) 
would  increase  the  difficulty  level  significantly  (see 
Francis  and  White  1974).  So,  although  a larger  number 
of  sources  will  increase  the  likelihood  of  obtaining  a 
final  solution  which  can  be  asserted  to  be  globally  optimal, 
it  is  likely  the  problem  will  be  of  such  a large  scale 
as  to  be  essentially  unsolvable  by  branch  and  bound. 
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Recognizing  the  difficulty  of  handling  such  large 
scale  problems,  Cooper  has  over  the  years  developed  a number 
of  heuristics  for  the  Euclidean  norm  problem  in  order  to 
get  a local  minimum.  Kuenne  and  Soland  (1971,  1972) 
coded  one  and  found  it  extremely  efficient.  A large  number 
of  local  minima  can  be  generated  rather  inexpensively. 

Since  the  single  facility  algorithm  for  the  spherical 
case  can  only  guarantee  a local  minimum  in  general,  it 
would  seem  practical  to  use  the  heuristic  for  the  sphere 
problem. 

The  basic  steps  of  the  Cooper  heuristic  used  by 
Kuenne  and  Soland  are: 

1.  Select  N initial  servicing  facility  locations. 

2.  Assign  the  M demand  points  to  the  closest 

source . 

3.  Solve  N single  facility  location  problems 
(results  in  a global  minimum  on  the  plane  for  each  sub- 
problem;  but  a local  or  global  minimum  on  the  sphere 
depending  on  whether  the  N subsets  are  each  contained  in  a 
disk  of  diameter  rr r / 2 ) . 

4.  Go  to  step  2. 

The  iterative  process  continues  until  some  convergence 
criterion  is  achieved. 

. For  practical  problems,  examination  of  a large 
number  of  local  minima  may  suffice.  Cooper  (1963)  found, 
as  did  Kuenne  and  Soland,  that  the  objective  function  is 
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fairly  fiat  and  hence  not  very  sensitive  to  changes  in 
location  around  the  global  optimum.  This  property  would 
appear  to  carry  forward  to  the  sphere.  In  fact,  one  could 
reasonably  expect  the  objective  function  to  be  flatter  on 
the  sphere.  On  the  sphere  the  single  facility  problem  can 
even  have  a constant  objective  function  value  for  ail 
points  in  the  domain.  This  occurs  when  six  equally 
weighted  points  are  placed  at  the  intersection  of  three 
orthogonal  great  circles.  Also,  the  upper  bound  of  ttR  on 
distances  between  points  limits  the  total  cost  at  any 
location.  This  is  not  so  for  the  planar  problem  since  the 
cost  can  increase  without  bound  by  selecting  progressively 
"worse"  locations  for  the  servicing  facility. 

It  is  important  to  note  that  use  of  the  heuristic 
can  also  help  in  the  branch  and  bound  approach  by  cheaply 
generating  a good  initial  feasible  solution  and  lower  bound. 

I I I . 1 4 . 5 The  Discrete  Multifacility  Problem 

In  addition  to  the  formulation  of  the  location- 
allocation  problem  given  by  (3.14.1),  other  formulations  have 
been  studied.  One  variation  involves  a discrete  solution 
space  and  is  known  as  the  plant  or  warehouse  location 
problem.  A treatment  of  this  problem  is  given  by  Francis 
and  White  (1974).  It  has  been  studied  by  ReVelle  and 
Swain  (1970),  Curry  and  Skeith  (1969),  and  Shannon  and 
Ignizio  (1970) . 
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It  is  of  interest  to  note  that  a discrete  location/ 


covering  problem  that  closely  parallels  the  objective  of 
the  foregoing  can  easily  be  handled  for  the  spherical 
case  by  using  existing  algorithms. 

An  efficient  heuristic  due  to  Shannon  and  Ignizio 
(1970)  will  provide  "good"  solutions  for  the  problem  on 
the  sphere.  The  fundamental  data  base  for  the  algorithm 
is  a fixed  matrix  [a.^]  of  distances/costs  derived  from 
interaction  between  discrete  servicing  facility  locations 
X.  and  demand  points  P^.  The  matrix,  of  course,  would 
incorporate  great  circle  distances,  thereby  making  it 
applicable  to  the  sphere  problem. 

Letting  W^p(X^,P^)  = a^  , the  objective  is: 

M 

Minimize  Z = E min  a. . 

i=l  LJ 

j £ e ( X) 

where:  0(X)  = {j|X.=l}^4> 

M 

E X.  < K 

j=l  J 

Xjtfo.i}  y. 

W j ^ is  the  weight  between  servicing  facility  X^ 
and  demand  point  P^ 


III. 15  Summary 

In  this  chapter  two  new  heuristic  algorithms  have 
been  presented  for  solution  of  the  single  facility  minisum 
location  problem  on  the  sphere.  In  addition,  several 
properties  of  the  minisum  sphere  problem  were  discussed, 
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and  extension  of  some  planar  results  to  the  sphere  were 
presented.  It  was  shown  that  the  search  for  a global 
optimum  may  be  restricted  to  the  spherically  convex  hull 
of  the  demand  points.  Application  of  the  solution  tech- 
niques to  solution  of  a special  class  of  multifacility 
problems  was  discussed. 


CHAPTER  IV 


MINIMAX  SINGLE  FACILITY  LOCATION  PROBLEMS 


IV . 1 Introduction 

In  this  chapter  minimization  of  the  maximum  dis- 
tance for  a single  facility  location  problem  on  a sphere 
in  the  great  circle  metric  is  addressed.  This  is  a natural 
extension  of  the  single  facility  minisum  problem  discussed 
in  Chapter  III  and  often  is  considered  to  be  a more  realis- 
tic objective.  Only  the  case  where  all  demand  point 
weights  are  equal  will  be  considered.  Two  similar  geo- 
metric approaches  to  the  minimax  problem  are  presented,  as 
are  a number  of  properties  of  the  problem. 

As  stated  in  Chapter  I,  in  an  norm  the  problem 
can  be  formulated  as : 


Minimize  Z 


Subject  to  |X-P.  L <_  Z 
1 P 


(A. 1.1) 


i=l M 

where:  X is  the  location  of  the  servicing  facility, 

P^  is  a demand  point  location,  i = 1 M,  and 

Z is  (geometrically)  the  radius  of  the  minimum 
covering  circle  or  sphere. 


90 


r 

i 

91 

As  observed  by  Elzinga  and  Hearn  (1972b)  the 
Euclidean  norm  (l^)  problem  can  be  solved  by  convex  program- 
ming techniques  for  a unique  global  optimum.  In  simplest 
form,  the  objective  is  to  find  the  smallest  circle  which 
covers  a finite  set  of  points  in  the  plane.  Capitalizing 
on  its  structure,  an  efficient  geometrical  technique  can 
also  be  used  to  solve  this  problem. 

IV . 2 Problem  Formulation 

On  the  unit  sphere  formulation  of  the  problem  is 
similar,  with  the  exception  that  the  great  circle  metric 
is  used,  and  a spherical  constraint  is  added.  Given  a 
finite  number  of  points  on  the  surface  of  the  sphere,  it 
is  desired  to  locate  a servicing  facility  on  the  sphere's 
surface  so  that  the  maximum  distance  between  this  point 
and  the  given  points  is  minimized.  That  is: 

Minimize  Z (4.2.1) 

Subject  to  Arccos  (a^x^-fb^x^+c^x^)  <■  Z i=l M 

2,2,2  -I 

Xj  + x.,  + x^=  1 

where:  (a^,b^,c^)  is  a demand  point,  i = 1 M 

(x^,x.,,x^)  is  a servicing  facility  location,  and 
Z is  (geometrically)  the  radius  of  the  minimum 
covering  spherical  cap. 


r 


I 

I 

I 


The  first  constraint  represents  the  distance  metric  and  the 
second  constraint  forces  the  solution  to  be  on  the  surface 
of  the  sphere.  The  second  constraint  unfortunately  creates 
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a non-convex  programming  problem,  as  will  be  seen  in 
Section  IV . 3 . 1 . 


It  is  easy  to  implicitly  eliminate  the  last  con- 
straint by  using  spherical  coordinates  and  formulate  the 

2 

problem  in  spherical  space  S : 

Minimize  Z (4.2.2) 

Subject  to  Arccos  (a^sin<t>cosA+b^sin4>sinA+c^cosif>) 

£ Z i = 1 M 


0 < < TT 

-ir  < A < it 

where  (sin^cosA  ,sin4>sinA  ,cos<}>)  is  the  location  of  the 

servicing  facility  and  (a.,b.,c^)  and  Z are  as  before. 

Recognizing  the  one  to  one  correspondence  between 

chord  and  arc  lengths  of  length  < it  the  problem  may  also 

3 

be  formulated  as  one  in  Euclidean  E space  by  adding  the 
constraint  which  forces  the  solution  to  be  on  the  sphere. 
That  is,  the  objective  is  to  equivalently  minimize  the  maxi- 
mum chord  length  (Euclidean  norm)  rather  than  minimize  the 
maximum  arc  length  (great  circle  metric) . It  would  be 
formulated  as : 

Minimize  Z (4.2.3) 

Subject  to  | (x]-ai)  2+(x^-bi)  ^+(x;i-ci)  2 1^  < Z 

2 , 2 , 2 , 

X1  + x2  + x3  = 1 


where : 


(a^,b^,c^)  is  a demand  point,  i = i,...,M 
(x^,x^,Xj)  is  the  servicing  facility  in  Euclidean 
three  space,  and 
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Z is  (geometrically)  the  radius  of  a sphere  with 
center  at  (x^.x^x^) 

IV. 3 Fundamental  Properties 

IV. 3.1  Non-convexity 

On  the  plane,  the  optimum  point  for  the  minimax 
problem  exists  and  is  unique.  This  is  intuitively  obvious 
and  is  shown  to  be  true  by  Blumenthal  and  Wahlin  (1941), 
among  others.  On  the  surface  of  a sphere,  however,  this 
uniqueness  cannot  be  guaranteed  for  the  general  case. 

As  observed  in  formulation  of  the  problem,  presence 
of  the  equality  constraint  results  in  a non-convex  program- 
ming problem.  This  is  most  vividly  demonstrated  by  example. 
Figure  4.1  depicts  the  objective  function  plotted  as  a 
function  of  colatitude  and  longitude  for  the  four  point 
problem  using  data  set  D5  of  Appendix  B.  Even  with  the 
small  number  of  points,  the  surface  is  full  of  well-defined 
ridges  and  valleys,  as  well  as  local  optima.  As  will  be 
pointed  out,  the  location  of  the  demand  points  on  the  sphere 
relative  to  one  another  is  all-important  in  determining 
how  easy  it  will  be  to  obtain  an  optimal  solution.  It  is 
intuitive  that  these  four  points  form  a global  subset  and 
are  not  containable  in  a hemisphere  (this  can  be  verified 
by  techniques  due  to  Blumenthal  (1936)).  Solution  of 
problems  having  demand  points  which  form  global  subsets 
cannot  be  handled  by  the  algorithms  presented  in  this 


A 
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chapter,  and  hence  remain  essentially  unsolved. 

Global  subsets  can  give  rise  to  alternate  optima-- 
which  cannot  occur  on  the  plane  due  to  the  convexity  of 
the  problem  in  Euclidean  space.  Perhaps  the  simplest 
example  is  the  case  when  six  demand  points  are  located  at 
the  intersection  of  three  orthogonal  great  circles.  In  such 
a case  there  are  eight  alternate  optima,  one  in  each  of  the 
octants  of  the  sphere  formed  by  the  intersecting  great 
circles  (see  Figure  4.2). 


* = visible  alternate 
optima,  with  four 
others  on  opposite 
hemisphere . 

Figure  4.2.  A special  case. 

A special  case  also  arises  when  the  set  of  demand 
points  includes  two  antipodal  points  and  the  set  is  not 
global.  In  such  a case  there  are  an  infinite  number  of 
alternative  optima.  They  consist  of  points  on  the  great 
circle  which  forms  the  equator  when  the  antipodal  points 
are  the  poles.  The  minimax  distance  is,  of  course, 
rr/2. 


A 
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IV. 3. 2 Comparison  with  En  Minimax  Problem 

The  n-dimensional  minimum  covering  sphere  problem 

has  been  investigated  by  Elzinga  and  Hearn  (1972b).  At 

first  glance,  it  might  appear  that  a solution  of  the  problem 
3 

in  E can  be  modified  to  solve  the  problem  on  the  sphere 
using  the  formulation  of  (4.2.3).  This  is  so  when  all 
demand  points  are  located  within  a hemisphere,  as  will  be 
seen,  but  it  cannot  be  used  when  the  demand  points  form  a 
global  subset.  For  example,  consider  again  the  six  points 
located  on  the  unit  sphere  S at  the  intersection  of  three 
orthogonal  great  circles.  The  minimum  covering  sphere 
(MCS)  is  the  sphere  S itself.  However,  the  covering  sphere 
determined  by  an  optimal  solution  to  the  minimax  problem  on 
the  sphere's  surface  is  much  larger  (see  Figure  4.3). 


Figure  4.3.  Minimum  covering  sphere  for  sphere 

problem. 

It  is  clear  that  for  the  problem  on  the  sphere's 


surface,  there  is  an  upper  bound  for  the  solution,  whereas 


there  is  none  on  the  plane.  This  bound  occurs  only  when 
all  points  on  the  sphere  are  demand  points.  Then  any  point 
on  the  sphere  is  optimal,  and  the  minimax  distance  is  n. 

In  the  problem  being  considered  there  are  a finite  number 
of  destinations  so  the  upper  bound  will  never  be  reached. 

Although  on  the  plane  two  points  determine  a unique 
circle  for  which  they  are  ends  of  a diameter,  on  a sphere's 
surface  the  points  can  determine  a great  circle,  a disk,  and 
a spherical  circle  whose  radius  is  greater  than  ^/2.  As 
will  be  seen,  when  working  with  global  subsets  these  lat- 
ter circles  are  the  ones  of  interest. 

IV. 4 Difficulties  with  Global  Subsets 

If  the  demand  points  cannot  be  contained  in  a 
hemisphere,  then  certainly  the  minimax  distance  is  greater 
than  tt / 2 , and  thus  the  covering  circle  is  a non-convex 
region  of  radius  greater  than  tt/2.  Many  difficulties  are 
encountered  in  attempting  to  solve  such  a problem,  some  of 
which  have  already  been  described. 

As  mentioned  in  Chapter  II,  a solution  was  offered 
for  the  sphere  problem  by  Sylvester  (1860)  , and  recognized 
more  recently  as  such  by  Blumenthal  and  Wahlin  (1941). 

The  solution  technique  was  actually  for  the  planar  problem, 
but  was  claimed  to  be  analogous  for  the  sphere. 

It  is  claimed  here  that  the  technique  outlined  by 
Sylvester  will  not  solve  the  problem  for  global  subsets. 
Sylvester  was  interested  in  the  problem  primarily  because 
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he  found  its  solution  useful  in  a paper  dealing  with  the 

n — 2 n — 2 — 2 

approximate  representation  of  ,/x  +y  and  vx  +y  +z  by 
linear  functions  of  x,  y and  z.  Early  in  the  paper, 
variables  are  restricted  to  be  positive.  It  appears  that 
this  may  have  led  to  the  error  in  concluding  that  the 
planar  technique  is  equally  applicable  to  the  sphere.  For 
example,  Sylvester  assumes  that  one  can  "...  reject  all 
those  points  that  are  contained  within  the  contour  formed 
by  the  arcs  joining  the  remaining  points,  so  that  the  case 
of  points  lying  at  the  angles  of  a convex  polygon  remain  to 
be  s tudied . . . . " 

Sylvester  did  not  consider  global  subsets  and  this 
is  where  the  approach  failed.  The  first  step  in  the  tech- 
nique is  to  find  a circle  (disk)  which  contains  all  the 
demand  points,  and  then  in  following  steps  progressively 
decrease  the  size  of  the  circle.  An  analogue  model  is 
offered  via  imagination  of  a rubber  band  which  maintains 
a circular  shape.  If  ".  . . sufficiently  stretched  over 

the  surface  of  the  sphere  to  contain  all  the  given  points 
(represented  by  minute  pins ' heads  given  upon  it) , this 
band  will  by  its  contraction  upon  the  surface  of  the  sphere 
. . . imitate  the  method  of  solution.  ...”  Clearly, 

this  procedure  will  not  work  for  global  subsets. 

An  efficient  technique  to  handle  the  minimax  objec- 
tive with  global  subsets  awaits  further  research.  Such  a 
problem  may  be  of  mathematical  interest,  but  it  will  only 
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arise  in  the  single  facility  case.  In  solving  a multifacil- 
ity  problem  there  is  no  need  to  consider  any  partition  of 
demand  points  which  forms  a global  subset.  Clearly,  when 
the  number  of  servicing  facilities  is  greater  than  or  equal 
to  two,  the  maximum  distance  between  any  demand  point  and 
its  servicing  facility  in  an  optimum  solution  is  bounded 
by  tt/2  merely  by  placing  the  servicing  facilities  at 
opposite  poles. 

Although  the  algorithms  presented  in  this  chapter 
deal  only  with  problems  for  which  the  demand  points  can  be 
contained  in  a hemisphere,  it  is  instructive  to  observe  the 
properties  of  the  single  facility  problem  with  demand  points 
forming  global  subsets  before  continuing. 

Once  it  is  known  that  the  demand  points  are  not 
containable  in  a hemisphere,  then  it  is  obvious  that  the 
minimax  distance  is  greater  than  tt/2.  Using  the  observation 
that  the  upper  bound  on  distances  on  the  sphere  is  tt  , it  is 
possible  to  convert  the  problem  to  an  equivalent  "maximin" 
problem. 

Consider  the  objective  of  finding  the  point  X^  which 
maximizes  the  minimum  distance  from  any  demand  point.  Sup- 
pose the  optimal  distance  is  b.  Consider  X^  , the  point 
antipodal  to  X^ . The  distance  to  any  demand  point  P^  from 
is  then  equal  to  ir-p (X^ , P^)  . Now  certainly  some  P^  is 
on  the  boundary  of  the  disk  D of  radius  b centered  at 
X^ . If  not,  D could  be  made  larger--a  contradiction.  So 
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for  this  , plX^.P^)  = b,  the  maximin  distance.  Now, 


P<W  = 


Max  Min  (p (Y , P . ) ) 

Y i 

Max[  Min(iT -p  (X,  P . ) ) ] , X antipodal  to  Y 
X i 1 

Max[  7i  -Min  (p  (X , P . ) ) J 
X i 1 

Max[  71  +Max(p  (X , P . ) ) ] 

X i x 


= 77  + Max! Maxp (X , P.  ) J 
X i 

= 77  - Min  Max(p(X,P.)) 

X i 1 

So  if  b,  the  maximin  distance  is  found,  one  has  the 
minimax  distance,  with  optimal  X^  being  antipodal  to  X^ . 

It  is  this  equivalence  which  illuminates  the 
difficulty  of  the  problem.  The  maximin  problem,  also  known 
in  the  literature  as  the  "noxious  facility"  problem  remains 
essentially  unsolved. 

The  problem  is  of  a combinatorial  nature  requiring 
consideration  of  the  different  disks  (less  than  tt  in  diameter) 
that  can  be  formed  using  the  demand  points.  The  objective 
is  to  find  the  largest  such  disk  which  does  not  contain 
any  demand  points  in  its  interior. 

As  an  example  of  the  foregoing,  consider  Figure  4. 4, 
which  depicts  two  problems  for  a degenerate  case  (on  a 
great  circle  arc).  In  Figure  A. 4(a)  the  solution  is  obvious 
and  can  readily  be  obtained  by  the  algorithm  outlined  in 
Sylvester's  paper.  In  Figure  4.4(b),  though,  there  is  no 
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disk  which  will  cover  all  the  demand  points,  and  Sylvester's 
approach  will  fail.  In  order  to  solve  it  intuitively,  it 
is  probably  easiest  to  follow  the  above  procedure.  That 
is,  find  the  largest  circle  (arc)  which  contains  no  demand 
points.  The  center  of  the  circle  (midpoint  of  the  arc), 

X*.  is  optimum  for  the  maximin  problem.  The  point  antipodal 
to  X*,  X*',  is  optimum  for  the  minimax  problem. 


Figure  4.4.  Global  subsets. 

I V . 5 Solving  the  Problem  on  a Hemisphere 
IV . 5 . 1 Introduction 

Recognizing  that  the  algorithm  of  Peirce  presented 
by  Sylvester  (1860)  will  solve  the  problem  when  all  points 
are  containable  within  a hemisphere,  the  question  naturally 
arises  as  to  the  value  of  any  other  approaches.  Peirce's 
approach  is  designed  for  solving  the  problem  by  hand  and 
is  not  amenable  to  programming  on  a computer.  Certainly 
if  the  only  concern  is  single  facility  problems,  though, 
this  procedure  is  quite  satisfactory.  However,  as  pointed 


102 


by  Hearn  (1971)  , high  speed  methods  would  be  desirable  for 
the  single  facility  problem  when  solving  multifacility 
problems.  Due  to  the  combinatorial  nature  of  the  multifacil- 
ity problem  (see  Section  IV. 6),  any  solution  procedure  is 
likely  to  involve  decomposition  into  single  facility  prob- 
lems, and  usually  the  resulting  number  of  such  problems  is 
quite  large. 

IV . 5 . 2 The  Euclidean  Norm  Technique 

Although  new  results  and  approaches  to  the  Euclidean 
norm  problem  are  continually  being  developed,  currently 
the  most  efficient  algorithm  for  solving  the  planar  problem 
is  a geometrical  approach  due  to  Elzinga  and  Hearn  (1972a). 
The  algorithm  begins  with  a circle  formed  by  any  two  points 
and  proceeds  to  monotonically  increase  the  size  of  the  circle 
until  all  points  are  contained  within  it.  Since  there  are 
only  a finite  number  ot  two  and  three  point  circles,  it  is 
a finite  algorithm. 

Elzinga  and  Hearn  (1972b)  also  provide  a solution 
procedure  for  the  problem  in  En.  It  involves  transforming 
the  Wolfe  dual  of  the  convex  programming  formulation  into 
a quadratic  programming  problem.  A finite  decomposition 
algorithm  based  on  the  Simplex  method  of  quadratic  program- 
ming is  developed. 

Encouraged  by  the  fact  that  Elzinga  and  Hearn's 


geometric  approach  is  currently  the  most  efficient  method 

2 

for  solving  the  problem  in  E , use  of  projective  geometry 
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was  made  to  capitalize  on  the  technique  in  order  to  solve 
the  problem  on  the  sphere.  The  only  additional  time 
required  to  solve  such  a problem  on  the  sphere  is  that  due 
to  the  transformation  to  the  plane  and  back  to  the  sphere. 

As  an  alternate  approach,  it  is  shown  that  the  solu- 
3 

tion  to  the  problem  in  E can  be  used,  via  a projection 
(normalization) , to  obtain  the  optimum  for  the  spherical 
problem.  In  practice,  however,  the  first  approach  would 
be  preferred  since  the  planar  algorithm  has  been  coded 
and  is  more  efficient  (Elzinga  and  Hearn  1972b)  . It  must 
be  kept  in  mind,  though,  that  both  approaches  require  that 
all  demand  points  be  containable  in  a hemisphere. 

I V . 5 . 3 The  Stereographic  Projection 

2 

The  geometric  algorithm  for  E can  be  effectively 
utilized  to  solve  the  sphere  problem  through  use  of  a stere- 
ographic projection.  This  projection,  most  likely  discovered 
by  Hipparchus  (circa  160-125  B.C.),  was  discarded  for 
centuries  when  the  world  was  decreed  flat. 

Kreysig  (1968)  credits  Lagrange  with  first  estab- 
lishing that  the  stereographic  projection  is  the  unique 
mapping  which  preserves  circles  from  the  sphere  to  the 
plane.  It  is  this  property  upon  which  one  can  capitalize, 
recognizing  the  geometric  basis  of  the  Elzinga-Hearn  algorithm. 

Cotter  (1966)  (cf.  Hilbert  1952)  provides  a simple 
proof  of  orthomorphism.  Referring  to  Figure  4.5  he  estab- 
lished the  similarity  of  AXab  and  AXAB.  He  then  observed 
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that  for  similar  triangles,  cones  of  rotation  are  similar 
and  the  projection  of  a circle  on  the  sphere  of  diameter 
AB  is  a circle  of  diameter  ab . 


Figure  4.5.  Stereographic  projection. 


Through  the  projection,  the  correspondence  is  one  to 
one  onto.  That  is,  to  each  circle  C on  the  plane,  there  cor- 
responds a unique  small  circle  C'  of  spherical  diameter  less 
than  tt  on  the  sphere. 

Now,  if  all  demand  points  are  contained  in  a hemi- 
sphere, and  the  projection  point  is  in  the  opposite 
hemisphere,  it  follows  via  the  projection  that  if  a point 
X is  contained  in  a circle  C on  the  plane,  its  correspond- 
ing inverse  image  point  X'  on  the  sphere  is  contained  in 
the  inverse  image  small  circle  C'  on  the  sphere.  Note 
that  this  is  NOT  true  if  all  points  cannot  be  contained 
in  a hemisphere  or  if  the  projection  point  is  not  in  a 
hemisphere  which  doesn't  contain  the  demand  points.  A 
suitable  projection  point  would  simply  be  a point  antipodal 
to  any  demand  point. 
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IV. 5.4  Stereographic  Projection/E  Technique 


Using  the  above  properties  and  observation  to  best 

2 

advantage,  a parallel  problem  can  be  constructed  in  E via 

the  stereographic  projection.  This  single  facility  unweighted 

minimax  problem  can  then  be  solved  in  tbe  plane  using  the 

2 

efficient  Elzinga-Hearn  algorithm.  In  E , the  optimum 
circle  C*  is  unique  (Rademacher  and  Toeplitz  1957) . Via 
the  projection,  C*  has  a unique  inverse  image  small  circle 
on  the  sphere. 

Requiring  all  demand  points  to  be  within  a hemisphere, 
and  appropriately  choosing  the  projection  point  in  the  oppo- 
site hemisphere,  note  that  the  spherical  small  circle  C*  con- 
tains all  the  demand  points  whenever  C*  contains  all  the 

P P 

images  of  the  demand  points.  Certainly  no  smaller  spherical 


circle  contains  all  the  demand  points.  C* , then,  is  the 
minimum  covering  circle  on  the  sphere.  Its  center  is  opti- 
mum for  the  spherical  minimax  problem  since  at  this  point 
the  maximum  distance  to  some  P.,  i=l,...,M,  is  obtained, 

l 

and  this  distance  can  be  no  smaller. 

Once  again,  it  is  important  to  note  the  necessity 
that  all  points  be  contained  in  a hemisphere,  otherwise 
the  inverse  transformation's  small  circle  corresponding 
to  the  optimum  circle  in  the  plane  will  not  contain  all 
demand  points  In  fact,  it  would  not  contain  any  demand 
points.  For  similar  reasons,  the  pole  for  projection  must 


be  exterior  to  the  spherically  convex  hull 
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demand  points.  That  is,  the  pole  must  be  in  a hemisphere 
which  contains  no  demand  points. 

3 

IV . 5 . 5 Normalization/E  Technique 

In  the  following  it  is  established  that  a solution 

technique  for  the  equally  weighted  single  facility  Euclidean 

3 

norm  minimax  problem  in  E can  be  advantageously  used  to 
solve  the  same  problem  on  the  unit  sphere  S if  all  demand 
points  are  on  an  open  hemisphere  of  S and  great  circle 
distances  apply. 

It  is  obvious  that  if  any  two  points  are  antipodal, 
or  if  the  cannot  be  contained  in  a hemisphere,  the  mini- 
mum covering  sphere  (MCS)  will  be  the  unit  sphere  S,  and 
the  optimum  point  would  be  the  center  of  S.  Such  a situ- 
ation precludes  solution  by  the  technique  which  follows. 

As  will  be  seen,  a crucial  part  of  the  procedure  is  to  pro- 
ject from  the  center  of  S through  the  center  of  the  MCS,  X* 
to  get  a unique  point  on  the  surface  of  S (normalization 
of  X*) . This  certainly  cannot  be  done  if  the  center  of 
S and  the  MCS  coincide.  A rudimentary  method  of  determining 
whether  all  points  are  in  a hemisphere  was  given  in 
Section  III . 4 . 3 . 

Prior  to  establishing  the  theorem  it  is  useful  to 
observe  that  on  a hemisphere,  the  correspondence  between 
the  small  great  circle  arc  and  the  chord  determined  by 
two  distinct  points  via  the  projection  of  points  of  the 
chord  onto  the  arc  from  the  sphere  center  is  one  to  one  onto. 
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Referring  to  Figure  4.6,  on  an  open  hemisphere 
it  is  easily  seen  that  as  X increases  monotonically  so  does 
A;  the  same  relationship  holds  between  x and  a.  It  also 
follows  that  a = R Arcsin(x/R) . 


Figure  4.6.  Arc  and  chord  relationship. 

Theorem  4. 1 Given  M equally  weighted  demand  points  on 

an  open  hemisphere  of  S,  suppose  the  Minimum  Covering  Sphere 
3 * 

in  E is  obtained.  Let  Xs  be  the  unique  projection  from 

the  center  of  the  sphere  S to  the  surface  through  X” , the 

*>v 

center  of  the  MCS . Considering  the  points  P^,  Xg  is  the 
optimum  for  the  minimax  distance  problem  for  the  great 
circle  metric  on  the  surface  of  S. 


Proof : Certainly  one  of  the  demand  points  will  be  on  the 

boundary  of  the  E MCS.  Otherwise  the  MCS  could  be  made 

smaller,  a contradiction.  Suppose  P ^ is  such  a point 

(see  Figure  4./).  Let  | X*  - P . | c be  the  optimum  minimax 

o j x o 

3 ^ 

distance  in  E , where  X*  is  the  center  of  the  minimum 

o 

covering  sphere.  Consider  the  chord  from  P,.  through  X* 

of  length  2 | X*  - P.|  = X.  Corresponding  to  this  chord 

o J *2 

is  a unique  small  great  circle  arc  of  measure  A.  Letting 
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X*  be  the  projection  of  X*  to  the  surface  of  the  sphere  it 

is  seen  that  p(X*,P.)  = \k. 

s J 

The  proof  that  X*  is  optimum  for  the  great  circle 
metric  on  S is  by  contradiction.  Let  Xg  be  a point  on  S 
distinct  from  X * such  that  its  maximum  distance  from  any 
other  point  is  p(X',P,),  where  P,  is  one  of  the  M demand 
points.  Further,  suppose  that  p(X',P,)  < p ( X* , P . ) . 

S K S J 

Now,  as  observed,  corresponding  to  p(X',P,  ) is 

length  | X^  - P^ | ^ via  the  inverse  correspondence.  It 

follows  that  | X ' - P,  | „ is  the  maximum  length  from  X'  to 
O K *2  o 

any  of  the  points  P^.  For,  if  not,  there  exists  P ^ sucn 


that  lx;  - Pli4> 


X' 

o 


P^  | ^ . Due  to  properties  of  the 


Arcsin  function  this  implies  that  p(X',Pj  > p(X',P,  ). 

s 1 s k 

This  cannot  be,  since  by  assumption  Pk  is  further  from  X'. 

P. 


In  the  same  manner,  the  assumption  that  p(X',P,)  < 

S K 


p(Xg,Pj)  implies  that  | X°  - Pk 


X*  - P . I „ . 
o J 1 


This  results  in  a conclusion  that  IX' 

' o 


pkU2  is  a 


maximum  distance  strictly  less  than  the  minimax  distance, 
a contradiction. 
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IV. 6 Application  to  Multifacility  Problem 
IV  6 . 1 Introduction 

As  mentioned  earlier,  the  primary  incentive  for 
developing  an  efficient  algorithm  for  the  single  facility 
problem  is  its  possible  usefulness  in  solving  a multifacil- 
ity problem.  In  this  section  a particular  type  of  multi- 
facility problem  is  considered  in  which  there  is  no  interac 

tion  between  facilities  to  be  located,  and  all  demand  point 

2 

weights  are  equal.  The  problem  in  E“~  was  briefly  discussed 
by  Hearn  (1971)  and  also  studied  by  Loginov  (1969)  for  the 
En  case.  The  mathematical  formulation  is: 

Minimize  V 

subj ect  to  max  {Z..*|X.-P.|»)  £ V i = 1 , . . . ,M 

P . , X . J1  J 1 2 j = 1 N 

i ’ J J 

N 

E Z . > 1 j = 1 N 

j = l J 

1 if  facility  X.  interacts  with 
Z.^  = { demand  poiAt  P. 

0 otherwise  1 * * * V 

where:  X.  is  the  location  of  a servicing  facility,  j = 1, . 

is  a demand  point  location,  i = 1,...,M 

V is  (geometrically)  the  radius  of  a disk  centered 
at  an  Xj , and 

M and  N are  the  number  of  demand  points  and  servic- 
ing facilities,  respectively. 

The  problem  can  be  stated  as:  Given  a finite  set 

A of  points  in  Euclidean  n-space,  and  an  integer  N, 


liO 


1 


determine  N n-dimensional  spheres  such  that  each  point  of 

A is  contained  in  at  least  one  of  the  spheres  and  the 

maximum  sphere  has  minimum  diameter.  It  can  be  considered 

as  a minimax  counterpart  to  the  multifacility  minisum 

location-al locat ion  problem  of  Section  III. 14. 

Loginov's  solution  procedure  will  solve  the  problem 
3 

for  small  N and  M in  E'  . Hearn  (1971)  presents  an  implicit 

enumeration  technique  which  will  solve  the  problem  in 
2 

E . Neither  technique  is  very  efficient,  but  at  least  a 
procedure  is  available. 

The  multifacility  minimax  problem  presents  the  same 
difficulties  encountered  in  the  multifacility  minisum  prob- 
lem. Due  to  its  combinator ial  nature,  current  solution  pro- 
cedures cannot  handle  more  than  3 or  4 servicing  facilities 
and  20  or  more  demand  points.  The  number  of  partitions  for 
any  given  M and  N is  S(M,N),  the  Stirling  number  of  the 

second  kind.  In  the  particular  case  of  N=2  , S(M,2)  = 

M-  1 

2 -1,  which  quickly  increases  for  moderate  size  M. 

The  intent  of  this  section  is  to  show  that  Loginov's 
3 

procedure  for  E"  , or  any  more  efficient  technique  that  may 
be  developed  in  the  future,  can  be  used  to  solve  the 
counterpart  problem  on  the  sphere. 

IV. 6. 2 Problem  Formulation 


For  the  problem  on  the  sphere,  the  formulation  is 
quite  similar,  except  that  great  circle  arcs  are  used,  and 
the  servicing  facilities  are  restricted  to  the  surface  of 


the  unit  sphere : 


Minimize  v 


Subject  to  Max{ p (X . , P . ) • Z . . } £ v i = 1 M 

J 1 J1  j = 1- N 


E Z . . > 1 

j=l  J1" 

|xih,  - 1 


i - 1 M 


j = 1, 


Z. . c { 0 , 1} 
J1 


V i.j 


Recognizing  that  as  long  as  all  demand  points  are 
contained  within  a hemisphere,  it  follows  that  the  maximum 
great  circle  arc  length  corresponds  to  the  length  of  the 
maximum  chord,  which  is  £ 2 for  the  unit  sphere.  Thus,  an 
equivalent  formulation  is  to: 


Minimize  v 


Subject  to  Max{ | X. -P. | 0 »Z..}  £ v 

J IX-,.  j i 


i = 1 M 

j = 1,  . . . ,N 


E Z. . > 1 

j-1  JL~ 

l*jU2  = 1 


1 = 1, 


i = l 


Zj  ie  { 0 , 1} 


V i.j 


IV. 6. 3 A Solution  Technique 


It  will  now  be  shown  that  any  procedure  which 
3 

solves  the  problem  in  E can  be  used  to  find  an  optimum 


solution  for  the  spherical  problem. 

Suppose  N £ 2.  As  mentioned  earlier,  the  optimum 
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minimax  distance  is  bounded  above  by  tt  / 2 . In  fact,  a mini- 
max distance  of  tt/2  would  occur  only  in  a pathological 
case  where  there  are  an  infinite  number  of  demand  facil- 
ities (a  great  circle  arc,  for  example).  In  this  instance 
there  are,  likewise,  an  infinite  number  of  alternate 
optima.  For  example,  any  two  antipodal  points  located 
on  the  great  circle  will  be  optimum  points. 

For  practical  purposes,  then,  it  can  be  assumed 
that  the  minimax  distance  is  strictly  less  than  ^/2.  This 
is  the  case  when  there  are  a discrete  number,  M,  of  demand 
points  and  N ^ 2. 

Suppose  all  points  P^,  i = 1,...,M  are  located  on  a 

3 

sphere  S.  S is  embedded  in  E . Capitalizing  on  this  fact, 
consider  the  problem  as  unconstrained  and  employing  an 
existing  algorithm  for  the  counterpart  multifacility  minimax 

3 

problem  in  E . 

3 

Solving  the  problem  in  E , the  solution  consists  of 
N spheres,  the  centers  of  which  are  servicing  facilities 
optimally  located  so  that  the  maximum  distance  to  any  demand 
point  is  minimized  and  all  P^  are  serviced  by  at  least  one 
facility.  The  union  of  the  N spheres,  then,  contains  all  the 
demand  points.  It  is  obvious  that  the  solution  need  not  be 
unique,  i.e.,  sometimes  it  is  possible  to  find  alternate 
systems  of  N spheres  covering  all  demand  points  such  that 
the  diameter  of  the  maximum  sphere  is  the  same.  In  any 
case,  of  course,  the  minimax  distance  is  unique. 


d 


| 
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Given  an  optimal  set  of  N spheres,  project  the 
centers  of  these  spheres  to  the  surface  of  S.  Since  no 
point  is  greater  than  distance  2 (diameter  of  the  unit 
sphere)  from  any  other  point,  the  N spheres  are  each 
certainly  no  larger  than  S.  Likewise,  as  mentioned  pre- 
viously, the  maximum  great  circle  distance  from  a projected 
point  to  any  demand  point  contained  within  the  respective 
sphere  is  less  than  u/2. 

Through  a similar  argument  as  in  the  single  facil- 
ity minimax  problem  for  points  on  a hemisphere  it  is  seen 
that  monotonicity  is  preserved  in  the  projection.  That 
is,  the  largest  sphere  generates  the  largest  great  circle 
distance.  So,  using  the  projected  points  as  the  N minimax 
locations  on  the  sphere,  the  largest  great  circle  distance 
is  found  by  projecting  from  the  center  of  the  largest  sphere. 

This  great  circle  distance  is  the  smallest  possible. 
For,  if  there  were  a smaller  one,  through  the  inverse 
transformation  there  would  exist  a sphere  with  diameter 
less  than  the  minimax  diameter,  a contradiction. 
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COMPUTATIONAL  RESULTS  AND  CONCLUSIONS 


V . I Introduction 

This  chapter  accomplishes  three  objectives.  It 
provides  a fundamental  comparison  of  the  two  algorithms  in 
Chapter  III  with  regard  to  their  ability  to  reach  a global 
optimum  (robustness) , provides  insight  to  the  properties  of 
the  single  facility  minisum  problem  on  the  spherical  sur- 
face, and  provides  evidence  of  the  error  which  arises  due 
to  the  use  of  a Euclidean  assumption  when  solving  large 
region  problems. 

A number  of  test  problems  have  been  solved  by  each 
algorithm,  varying  the  level  of  difficulty.  Two 
example  "large  region”  problems  from  the  literature  have 
also  been  solved.  In  each  problem,  verification  of  global 
and  local  optima  was  made  by  use  of  three  dimensional 
plotting  and  an  exhaustive  grid  search. 

Comparisons  are  made  between  solutions  found  via  a 
Euclidean  assumption  and  those  found  using  a great  circle 
metric.  The  Hyperboloid  Approximation  Procedure  (HAP) 
(Eyster,  et  al.  1973)  is  used  to  solve  the  problems  for 
the  Euclidean  assumption.  It  is  seen  that  even  in  the 
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cases  where  the  optimum  location  is  not  very  different, 
the  difference  in  objective  function  values  can  be  signif- 
icant. Other  conclusions  are  also  drawn  based  upon  compu- 
tational experience. 

The  two  single  facility  algorithms  of  Chapter  III 
were  coded  in  Fortran  IV,  and  all  problems  were  solved  on 
the  IBM  370/158J  system  using  double  precision  arithmetic. 
No  attempt  was  made  to  determine  sensitivity  of  parameter 
selection.  For  example,  judicious  selection  of  the  step 
sizes  could  significantly  affect  the  ability  of  the  Cyclic 
Meridian  Search  (CMS)  algorithm  to  reach  global  optima,  but 
confirmation  is  a subject  for  future  research.  Also, 
neither  algorithm  was  coded  with  the  intent  of  achieving 
high  efticiency.  The  primary  goal  was  to  reach  a global 
optimum.  For  more  efficient  codes  it  is  necessary  to 
incorporate  acceleration  procedures. 

The  parameters  provided  in  Table  5.1  were  used 
for  solving  the  example  problems  in  Appendix  B under  the 
great  circle  metric.  Refer  to  Chapter  III  for  descriptions 
and  definitions  of  each  parameter. 
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Planar  Projection 

Cyclic  Meridian  Search 

BOUND  = tt/6 

Al  = 10"1 

a =2 

EP1  = 10"6 

AINC  = 10-1 

ITBOUN  = 40 

EP2  = 10-6 

AMIN  = 10“6 

LRANBO  = 4 

N1  = 50 

e = 10"9 

S =1 

— 

EP1  = 10-6 
EP2  = 10“6 
N =5 

Table  5.1.  Parameter  values. 

V . 2 Some  Example  Problems 

The  first  example  is  a minisum  problem  using  the 


demand  points  in  Data  Set  D6  of  Appendix  B.  It  consists  of 
13  equally  weighted  points  located  on  a great  circle  arc.. 

It  is  informative  since  its  simple  nature  makes  it  easy  to 
see  the  effect  that  different  starting  solutions  have  on 
the  iterative  procedures  of  both  algorithms.  Figure  5.1 
depicts  a graph  of  the  objective  function  values  as  a func- 
tion of  longitude.  The  latitude  is  0°  for  all  demand 
points.  The  demand  points  are  indicated  in  the  figure, 
as  are  the  locations  of  all  starting  solutions  except  S5. 
Note  that  a convergent  algorithm  would  be  expected  to 
achieve  the  global  optimum  (0,20)  for  any  starting  solution 
between  -60°  and  56°  longitude.  Both  algorithms  behave 
properly,  with  a surprising  result  arising  from  use  of  CMS. 
Using  starting  solution  S2,  one  would  expect  to  get  stuck 
at  a local  minimum,  yet  due  to  the  internal  features  of 
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the  algorithm  (discussed  in  Chapter  III)  the  search  escapes 
the  local  minimum's  region  of  attractiveness,  jumps  across 
another  one,  and  stops  at  the  global  optimum.  This  feature 
seems  to  occur  rather  frequently. 

Also,  note  that  when  starting  solution  S5  is  used, 
the  Planar  Projection  Algorithm  (PPA)  reaches  the  global 
optimum  while  CMS  stops  at  a local  optimum.  Curiously,  in 
all  example  problems  attempted  this  was  the  only  instance 
in  which  PPA  achieved  a global  optimum  when  CMS  did  not. 

It  is  enough,  however,  to  encourage  one  to  use  each  algorithm 
as  a check  on  the  other. 

In  Table  5.2  the  results  of  trying  a number  of  dif- 
ferent starting  solutions  on  Data  Set  D6  are  listed.  Note 
that  in  this  problem  a projected  centroid  starting  solution 
(S4)  results  in  a local  minimum,  regardless  of  the  technique 
used . 


Starting 

# 

Stops 

Obj  Fn 

Type 

# 

Stops 

Obj  Fn 

Type 

Point 

Iter 

at 

Value 

Soln 

Iter 

at 

Value 

Soln 

SI (0 , 20) 

1 

(0,20) 

17.52 

Global 

1 

(0,20) 

17.52 

Global 

S2 

(0,-155) 

3 

(0,-124) 

22.69 

Local 

2 

(0,20) 

17.52 

Global 

S3  (0,0) 

2 

(0,20) 

17.52 

Global 

1 

(0,20) 

17.52 

Global 

S4  (0, 
62.27) 
(Proj  . 
Cent . ) 

2 

(0,70) 

17.91 

Local 

1 

(0,70) 

17.91 

Local 

S5 

(90,90) 

9 

(0,70) 

1/.52 

Global 

6 

(0,70) 

17.91 

Local 

Table  5.2.  Results  for  Data  Set  D6 
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Now  consider  a three  equally  weighted  demand  point 
problem  using  Data  Set  Dl  of  Appendix  B.  This  problem 
tests  the  ability  of  the  algorithms  to  find  a global  optimum 
over  a definitely  non-convex  surface.  A three  dimensional 
function  plot  of  the  problem  is  found  in  Figure  3.1.  The 
global  optimum  is  at  P3,  (30,-60). 

The  algorithms  stopped  at  the  same  point,  as  can 
be  seen  in  Table  5.3,  for  all  but  one  of  the  starting  solu- 
tions. An  interesting  observation  is  that  when  one  starts 
at  the  antipodal  point  to  the  global  optimum,  P3',  both 
searches  still  end  up  at  the  global  optimum.  Concerning 
the  number  of  iterations,  in  general  CMS  is  more  efficient. 
Although  not  concerned  with  efficiency  in  terms  of  execu- 
tion times  in  this  research,  it  is  interesting  to  note 
that  CMS  takes,  on  the  average,  about  half  the  time  as  PPA. 
There  are  situations,  though,  when  CMS  is  slower  than  the 
other  approach,  particularly  when  the  search  begins  on  or 
near  a ridge. 

The  next  set  of  examples  consists  of  ten  small 
problems  for  which  a projected  centroid  was  used  as  a 
starting  solution.  As  can  be  seen  in  Table  5.4,  both 
algorithms  successfully  found  the  global  optimum  in  all  but 
one  case. 

This  group  of  problems  gave  rise  to  a number  of 
interesting  observations.  For  example,  Data  Set  D12  was 
handled  more  efficiently  by  PPA.  Data  Set  D9 , on  the 
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Data 

Set 

# 

Pts 

Type 

Soln 

No. 

Iterations 

Proj  . 
Centroid 

Obj 

Fn 

Value 

Optimum 

Soln 

Obj 

Fn 

Value 

PPA 

CMS 

D3 

7 

Global 

3 

4 

(25.2,27.6) 

6.46 

(25.8,21) 

6.15 

D7 

3 

Global 

3 

H 

(32.3,0) 

2.8 

(54.9,0) 

2.59 

D8 

3 

Global 

3 

fl 

(32.3,-60) 

2.8 

(54.9,-60) 

2.59 

D9 

3 

Global 

6 

a 

(63.6,75) 

3.64 

(80.6,75) 

DIO 

3 

Local 

7 

3 

(0,90) 

4.11 

(-5.3,90) 

4.11 

Dll 

3 

Global 

3 

5 

(27.4,71.5) 

2.78 

(50,80) 

2.6 

D12 

4 

Global 

2 

7 

(/0.8.37 .8) 

1.44 

(80,45) 

1.3 

D13 

4 

Global 

6 

3 

(49.3,110) 

4.41 

(44.8,50.3 

) 4 . 23 

D14 

4 

Global 

5 

3 

(45.9,88.3) 

3.79 

(52.1,100) 

3.76 

D15 

6 

Global 

1 

1 

(0,0) 

9.42 

(.97,0) 

9.42 

Table  5.4.  Results  for  ten  problems. 


other  hand  was  dispensed  with  in  one  iteration  by  CMS.  This 
is  because  in  the  first  iteration  the  search  occurred  along 
the  meridian  upon  which  the  optimum  solution  exists. 

In  the  large  majority  of  cases,  the  first  iteration 
makes  a major  move  toward  the  global  (or  local)  optimum,  and 
the  remaining  iterations  involve  small  steps  which  continu- 
ally, but  slowly,  improve  the  objective  function  value  until 
a stopping  criterion  is  satisfied.  In  order  to  illustrate 
this  property  of  the  algorithmic  search  pattern,  consider 
Data  Set  D7  as  solved  by  PPA: 
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Location 

Objective  Function 
Value 

Starting  Solution: 

(32.3,0) 

2.797699 

Iteration  1: 

(53.97,0) 

2.59099 

Iteration  2: 

(54.84,0) 

2.590403 

Iteration  3: 

(54.86,0) 

2.590402 

Certain  observations  about  some 

of  the  ten  Data 

Sets  are  of  interest  in  light  of  comments  and  results  in 
Chapter  III. 

1.  Data  Sets  D3  and  D12  are  containable  in  a disk 
of  diameter  tt/2,  and  can  be  considered  via  Conjecture  3.7.1. 

2.  Data  Set  D14  illustrates  Fagnano's  Result 
extended  to  the  sphere. 

3.  Data  Sets  D/  through  Dll  demonstrate  aspects 
of  Steiner's  Problem  on  the  sphere.  It  can  be  easily 
verified  that  the  optimum  solution  is  at  a point  which 
forms  angles  of  at  least  120°  with  the  three  vertices  of 
the  spherical  triangles  formed  by  the  demand  points.  The 
global  optimum  for  Data  Set  Dll  occurs  at  a vertex,  which 
has  an  angle  greater  than  120°.  There  are  two  alternate 
global  optima  for  Data  Set  DIO,  but  neither  algorithm 
reached  one  by  using  the  projected  centroid  as  a starting 
solution.  The  alternate  optima  occur  at  the  vertices 
(-30,20)  and  (-30,160).  The  Local  optimum,  at  which 
iterations  stopped,  is  the  interior  point  of  the  spherical 
triangle  which  forms  angles  of  precisely  120°  with  lines 


drawn  to  the  three  vertices. 
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4.  Data  Set  D15  has  a constant  objective  function 
value.  Every  point  on  the  sphere  is  optimal. 

Now  consider  a problem  using  Data  Set  D16  of 
Appendix  B.  This  six  unequally  weighted  demand  point 
problem  was  solved  by  both  algorithms  using  a variety  of 
starting  solutions  consisting  of  the  optimum  found  via  a 
Euclidean  assumption,  the  projected  centroid,  all  demand 
points,  the  antipodal  point  to  the  projected  centroid  and 
antipodal  points  to  all  demand  points.  By  referring  to 
Figure  5.2,  one  can  see  intuitively  the  expected  "goodness” 
or  "badness"  of  the  respective  starting  solutions  relative 
to  the  location  of  the  global  optimum. 

The  computational  results  are  exhibited  in  Table 
5.5.  For  the  various  starting  solutions  used,  CMS  per- 
forms better  overall  than  PPA  in  reaching  the  global  optimum. 
Also,  in  execution  time  CMS  achieved  the  solution  twice  as 
fast,  on  the  average,  as  PPA. 

Of  particular  interest  are  the  results  when  starting 
solutions  S6  and  S8  are  used.  Although  they  are  local 
minima,  CMS  successfully  escapes  their  grasp  and  reaches 
the  global  optimum. 

Note  that  in  this  problem  the  projected  centroid  is 
a good  starting  solution,  and  is  not  too  far  from  the  optimal 
solution.  In  fact,  it  is  significantly  better  than  what 
would  have  been  achieved  if  the  problem  had  for  some  reason 
been  solved  with  a Euclidean  assumption.  Table  5.6 
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presents  the  results  for  the  projected  centroid,  the  optimum 

for  the  sphere,  and  the  optimum  found  by  using  a Cartesian 

2 

coordinate  system  in  E . 


Location 
(Latitude , 
Longitude) 

Objective  Function 
Value 

Great  Circle 

Op  t imum 
for  Sphere 

(25,-115) 

20.57 

37.486 

Optimum 
for  E2 

(9.44,21.02) 

22.195 

25.214 

Proj  ected 
Centroid 

(42.4,-91) 

20.90 

35.434 

Table  5.6.  Sphere  vs.  plane.  Data  Set  D16. 

Using  the  Euclidean  assumption  for  this  particular 


problem  will  result  in  an  objective  function  value  error  of 
about  18%  and  a location  error  of  over  10 , 000  miles ! 

The  next  example  required  the  algorithms  to  solve 
a problem  in  the  polar  regions.  Data  Set  D17  of  Appendix  B 
includes  the  coordinates , in  latitude  and  longitude  to  the 
nearest  degree,  of  nine  experimental  stations  and  other 
sites  in  Antarctica.  Using  the  minisum  criterion  and 
assuming  equal  weights,  an  optimal  solution  was  found  by 
both  algorithms.  The  projected  centroid  was  used  as  a 
starting  solution.  PPA  used  three  iterations  while  CMS 
used  four  iterations  to  reach  the  optimum  at  (-87.79,89.78). 

With  the  foregoing  observations  in  hand,  two  problems 
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from  the  literature  were  examined.  Both  involved  solution 
of  a single  facility  minisum  problem  by  using  a Euclidean 
assumption  and  letting  latitude  and  longitude  represent  the 
Cartesian  coordinates. 

First,  consider  an  example  problem  used  by  Kuhn  and 
Kuenne  (1962).  The  data  (Data  Set  D18  of  Appendix  B)  con- 
sists of  24  cities  of  the  Ukrainian  Soviet  Socialist 
Republic.  The  weights  are  equal  to  the  proportion  each  had 
of  the  population  of  the  100  largest  Russian  cities  (as 
of  1961) . The  coordinates  are  degrees  of  north  latitude 
and  east  longitude  correct  to  the  nearest  full  degree.  The 
latitude  of  Kharkov  is  in  error  by  about  nine  degrees , 
but  since  Kuhn  and  Kuenne  evidently  used  the  published 
coordinates  to  solve  the  problem,  they  were  also  used  for 
this  example. 

Before  examining  the  results,  it  is  important  to 
point  out  that  the  points  could  easily  be  fit  into  an  800 
by  1000  mile  rectangle,  and  certainly  into  a disk  of  diameter 
less  than  ttR/2,  where  R is  the  radius  of  the  earth.  By 
Conjecture  3.7.1,  one  can  guarantee  a global  optimum.  Both 
algorithms  achieved  the  global  optimum  using  the  projected 
centroid  as  a starting  solution. 

In  this  particular  problem  there  is  no  significant 
effect  on  location  of  the  optimal  solution  relative  to 
using  either  the  great  circle  metric  or  Euclidean  norm. 
However,  there  is  a definite  effect  on  the  value  of  the 
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objective  function  (see  Table  5.7).  This  would  be  of  sig- 
nificant interest  in  the  area  of  logistics,  and  especially 

(budgeting,  in  any  attempts  to  estimate  transportation 
requirements/ costs . 


Location 
(Latitude , 
Longitude) 

Objective  Function 
Value 

Great  Circle 

Euclidean 

Optimum 
for  Sphere 

(47.795,35.2402) 

.01293768 

.01669464 

2 

Optimum  E 
(HAP) 

(47.707,35.104) 

.01295154 

.01663096 

2 

Optimum  E 
(Kuhn/Kuenne) 

(47.6,35.32) 

.01296313 

.0166528 

Table  5.7.  Sphere  vs.  plane,  Data  Set  D18. 

It  is  easy  to  verify  that  the  error  in  optimal 


location  is  about  9.26  miles,  while  the  objective  function 
value  error  is  off  by  about  227o. 

The  last  example  is  due  to  Chapelle  (1969)  . He 


considered  the  49  continental  states  plus  the  District  of 
Columbia  and  their  corresponding  output  in  first  class 
letter  mail  in  the  year  1965.  The  output  volume  in  pounds 
was  assumed  to  originate  from  the  capital  of  each  state. 
Latitude,  longitude  and  weight  for  each  point  location  are 
found  in  Data  Set  D19  of  Appendix  B.  Chapelle's  objective 
was  to  optimally  locate  a central  facility  for  the  postal 
system,  assuming  plane  geometry. 
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Although  the  region  is  larger  than  in  the  previous 
example,  it  still  can  be  contained  in  a disk  of  diameter 
less  than  ttR/2.  Both  algorithms  of  Chapter  III  reached 
the  global  solution,  using  the  projected  centroid  as  a 
starting  solution,  taking  seven  iterations  of  CMS  and  three 
iterations  of  PPA  to  meet  stopping  criteria.  It  took  657, 
less  time  to  solve  the  problem  by  CMS. 

Referring  to  Table  5.8,  one  can  compare  the  results 
obtained  by  either  a Euclidean  norm  or  great  circle  metric 
approach.  The  error  in  location  is  only  30  miles,  yet 
the  error  in  objective  function  value  is  over  207,. 


Location 
(Latitude , 
Longitude) 

Objective  Function 
Value 

Great  Circle 

Euclidean 

Optimum 
for  Sphere 

(40,83) 

780,065 

987,413 

Optimum 
for  Plane 

(39.65,-83.45) 

781,635 

986,878 

Table  5.8.  Sphere  vs.  plane,  Data  Set  D19 . 

V . 3 . Conclus ions 

Based  upon  the  experience  obtained  in  solving  a 


number  of  single  facility  minisum  problems  using  both 
algorithms,  the  following  conclusions  can  be  drawn: 

1.  Comparing  the  two  algorithms  in  their  rudi- 


mentary form,  the  Cyclic  Meridian  Search  (CMS)  algorithm 
is,  in  general,  faster  than  the  Planar  Projection 
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Algorithm  (PPA) . There  is  always  the  possibility  of 
getting  stuck  on  a ridge  and  zigzagging,  but  the  integral 
use  of  random  search  directions  helps  avoid  these  diffi- 
culties . 

2.  PPA  is  not  as  likely  to  reach  the  global  opti- 
mum as  CMS.  Other  than  the  observation  that  CMS  can  "jump 
out"  of  local  minima  regions  and  accelerate  in  a search 
direction,  there  is  no  apparent  a priori  reason  for  its 
success . 

3.  No  case  of  divergence  was  observed  for  either 
algorithm  in  any  problem  attempted.  In  all  cases  at  least 
a local  minimum  was  obtained. 

4.  The  projected  centroid  is  a good  starting  solu- 
tion, but  will  not  necessarily  result  in  finding  the  global 
optimum. 

5.  A Euclidean  assumption  for  even  moderately 
large  regions  (800  by  1000  miles)  can  lead  to  significant 
error  in  objective  function  values.  This  is  because  a 
planar  assumption  overestimates  actual  distance.  Also, 
significant  error  in  location  of  the  optimal  facility  can 
occur  in  problems  covering  regions  as  large  as  a hemisphere. 


CHAPTER  VI 


SUMMARY  AND  RECOMMENDATIONS  FOR  FUTURE  RESEARCH 


VI . 1 Introduction 

In  this  chapter  a summary  of  research  performed  on 
the  single  facility  location  problem  for  large  regions  is 
presented.  Recommendations  for  future  research  are  included 
to  assist  in  identifying  topics  which  would  extend  the  cur- 
rent research  effort. 


VI . 2 Summary 

This  research  presents  an  extension  of  the  gen- 
eralized Weber  problem  and  of  the  related  problem  of  mini- 
mization of  maximum  distance.  The  primary  objectives  were 
to  investigate  the  effect  of  a great  circle  metric  on  the 
problems  of  interest,  and  to  develop  solution  procedures 
to  solve  them.  The  results  will  permit  greater  realism 
and  accuracy  in  solving  large  region  location  problems 
for  which  a Euclidean  (planar)  assumption  is  inappropriate. 

In  Chapter  I the  necessity  of  considering  location 
problems  under  a great  circle  metric  is  discussed,  along 
with  possible  applications  of  the  research.  Chapter  II 
then  presents  a literature  survey  of  previous  research 
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concerning  large  region  location  problems. 

In  Chapter  III  the  single  facility  generalized 
Weber  problem  is  studied,  where  the  goal  is  to  locate  a 
new  facility  on  the  sphere  relative  to  several  weighted 
demand  facilities  such  that  the  total  cost  of  transporta- 
tion is  minimized.  Cost  is  considered  to  be  strictly  a 
function  of  great  circle  distance  and  demand  point  weights. 

Different  formulations  and  properties  of  the  problem 
were  discussed,  along  with  two  new  iterative  algorithms 
for  solving  it.  Both  algorithms  were  programmed  and  optimal 
solutions  were  obtained  for  several  problems.  Application 
of  the  techniques  to  solving  a multifacility  location- 
allocation  problem  is  also  discussed. 

In  Chapter  IV  the  objective  of  locating  a facility 

to  minimize  the  maximum  distance  to  any  demand  point  under 

a great  circle  metric  was  considered.  It  is  assumed  that 

the  number  of  demand  facilities  is  finite  and  their  weights 

are  all  equal.  Two  approaches  to  the  problem  are  presented, 

both  of  which  capitalize  on  existing  algorithms  to  solve 

2 3 

the  Euclidean  norm  problem  in  E or  E Properties  of  the 
problem  and  application  of  the  single  facility  techniques 
to  solution  of  a certain  multifacility  minimax  problem  are 
presented . 

Chapter  V presents  results  of  applying  the  algorithms 


in  Chapter  III  to  a number  of  example  problems.  Comparative 
effectiveness  in  reaching  the  global  optimum  is  examined. 
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VI. 3 Recommended  Future  Research 

A number  of  areas  for  further  study  are  clearly 
available.  They  are: 

1.  Perform  sensitivity  analyses  for  the  two 
algorithms  in  Chapter  III  and  compare  their  relative 
efficiency . 

2.  Due  to  the  nature  of  the  minisum  problem,  an 
efficient  global  algorithm  awaits  a breakthrough  in  non- 
convex  programming.  In  the  meantime,  efforts  toward  devel- 
oping good  bounds  using  a branch  and  bound  technique  would 
be  fruitful. 

3.  With  regard  to  escaping  local  optima  in  the 
minisum  problem,  investigate  development  of  existing  tech- 
niques to  capitalize  on  the  structure  of  the  problem. 

4.  Develop  an  efficient  technique  to  solve  the 
single  facility  minimax  problem  with  equally  weighted 
demand  points  when  the  demand  points  form  a global  subset. 
Based  on  comments  in  Chapter  IV,  success  in  this  endeavor 
awaits  a breakthrough  in  solution  of  the  "noxious  facility" 
location  problem. 

5.  Develop  a method  to  solve  the  weighted  single 
facility  and  multifacility  minimax  problem  in  the  great  circle 
metric.  it  does  not  appear  likely  that  the  geometric 
approach  would  be  feasible  for  such  problems. 

6.  Recognizing  that  in  practical  applications  one 
is  faced  with  the  fact  that  land  covers  only  a small  part 
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of  the  earth,  develop  methods  which  will  consider  con- 
straints as  to  location  of  the  servicing  facility  (besides 
the  one  restricting  it  to  the  sphere's  surface). 
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APPENDIX  A 


CONVEXITY  THEORY  IN  SPHERICAL  GEOMETRY 


A. I Preliminary  Definitions 

In  spherical  geometry  one  has  the  axiomatic  frame- 
work as  set  forth  by  Kay  (1969)  with  the  canonical  model 

2 

being  the  unit  sphere  in  Euclidean  3-space  given  by  |xj„  = 

2 

1.  In  this  model,  the  spherical  metric  p(x,y)  for  each 
x and  y in  S"^  is  the  length  of  the  minor  arc  of  the  great 
circle  joining  x and  y,  and  a,  the  least  upper  bound  for 
the  metric,  has  the  value  tt  . Further,  the  lines  of  the 
geometry  are  great  circles,  and  if  two  points  x and  y are 
not  antipodal  [p(x,y)  ^ tt  ] the  line  through  x and  y is 
unique;  in  this  case  the  line  is  denoted  xy . Each  line  can 
be  coordinatized  from  the  set  of  numbers  (A|-tt  < A < tt} 

(the  coordinate  set)  , with  any  point  as  origin  (zero  coordi- 
nate) and  such  that  if  the  coordinates  of  u and  v are  A and 
u,  then  p(u,v)  = | A - p | if  |A-p|  £ tt  , or  p(u,v)  = 2tt  - 
| A - u | it  | A - m | > 7i.  betweenness  can  be  defined  by 
asserting  that  z is  between  x and  y,  if  and  only  if  x,  y, 
and  z are  distinct  points  and  p(x,z)  + p(z,y)  = p(x,y) . 

If  z is  between  x and  y we  denote  the  relationship  by  (xzy) . 
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For  any  two  points  x and  y such  that  0 < P(x,y)  < v 
the  segment  joining  x and  y,  denoted  xy , is  the  set  of 
points  {z|z=x,  z=y,  or  (xzy)}. 

2 

Definition  A . 1 A set  C in  S is  said  to  be  convex  if 
for  each  two  points  x and  y in  C such  that  0 < p(x,y)  < it 
then  the  segment  xy C C.  The  convex  hull  of  a set  A in  S'" 
is  the  intersection  of  all  convex  sets  containing  A, 
denoted  conv  A. 

It  is  obvious  that  an  arbitrary  intersection  of 

convex  sets  is  convex;  hence  conv  A is  a convex  set  for  each 

subset  A of  S^.  A set  is  convex  iff  it  equals  its  own  convex 

hull.  Note  that  the  above  definition  of  convexity  allows 

the  following  examples  to  be  convex  sets:  A hemisphere,  a 

great  circle,  a lune , and  any  pair  of  antipodal  points. 

The  convex  hull  of  2 points  is  either  a segment  or  the 

two  points  themselves;  the  convex  hull  of  a closed  hemisphere 

2 

and  a point  not  on  it  is  the  entire  sphere  S . 

One  of  the  important  axioms  concerning  betweenness 

2 

on  S is  the  so-called  plane  separation  property:  Each 

2 2 

line  L in  S divides  the  points  of  S into  three  pairwise 
disjoint  convex  sets  L,  H^,  and  H.,  such  that  if  xr.H^  and 
there  exists  a point  zeL  such  that  (xzy) . The  sets 
H| , are  referred  to  as  the  ha  1 f-spaces  determined  by  L, 
and  two  points  which  lie  in  the  same  half-space  determined 
by  L are  said  to  lie  on  the  same  side  of  L.  Note  in 
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2 

passing  that  a half-space  is  a hemisphere  on  S with  its 
boundary  (great  circle)  deleted. 

A half- line  or  ray  is  in  a subset  L ' of  a line  L of 
the  form  {u[A]|X>_0}  where  u[A]  denotes  a coordinate  system 
for  L;  a half-line  always  consists  of  a (closed)  semi- 
great circle.  The  point  ufOl  = x is  called  the  origin  of 
the  ray,  and  if  y is  a point  on  L'  such  that  0 < p(x,y)  < n , 
the  ray  is  denoted  1$.  It  is  shown  by  Kay  (1969)  that  if 
zexy  with  0 < p(x,z)  < it,  then  xy  = "x£. 

An  angle  is  the  union  of  two  rays  having  the  same 
origin,  the  two  rays  being  called  its  s ides  and  the  common 
origin  its  vertex . If  xy  and  x£  are  the  sides  of  an  angle 
the  notation  -$.yxz  is  used.  The  measure  £yxz  of  the  angle 
£ yxz  is  the  Euclidean  measure  of  the  angle  formed  by  the 
tangents  to  the  sides  "xy  and  xz  at  x in  3-space.  Note 

2 

that  a straight  angle  (having  measure  if)  is  a line  on  S , 

and  that  every  other  angle  is  the  boundary  of  a region 
2 

on  S usually  called  a lune  or  spherical  section . 

Given  three  distinct  rays  xu , xv,  having  the 

same  origin  x,  we  say  that  "x^  lies  between  xu  and  ”xw , and 

we  write  (xu  x\^  xw)  , if  4 uxv  +-$.vxw  = ^(.uxw.  The  interior 

of  an  angle  -£yxz  is  the  set  of  all  points  lying  on  rays 

which  lie  between  the  sides  "xy  and  xz . (See  Kay  (1969) 

for  betweenness  properties  of  segments,  rays,  and  triangles.) 

2 

If  x,  v,  z are  3 distinct  noncollinear  points  of  S 
it  follows  that  0 < P(x,y)  < n , 0 < P(x,z)  < u , and 
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0 < P(,y,z)  < tt  . Thus,  the  segments  "xy , "yz , "zx  exist,  their 
union  defining  a set  called  a triangle , denoted  Axyz . If 
for  each  point  x and  line  L not  through  it  H(x,L)  denotes 
the  haif-plane  determined  by  L and  containing  x,  the 
interior  of  triangle  xyz  is  by  definition  the  set  H(x,yz)ft 
H(y  ,£?)r» H(z  ,1c^)  . It  can  be  shown  that  p is  an  interior 
point  of  triangle  xyz  iff  the  line  Icp  intersects  yz  at  a 
point  w such  that  (xpw)  and  (ywz).  and  similarly  with  lines 
y£  and  tp.  The  Crossbar  Principle  (Kay  1969)  then  shows 
that  (xy  3c^  x£)  . 

2 

Finally,  a circle  on  S is  the  set  of  all  points  at 
a fixed  distance  (radius)  from  a fixed  point  (center) . 

2 

Note  that  a circle  with  radius  it/ 2 is  also  a line  on  S ; 

the  center  of  any  other  circle  may  be  chosen  in  such  a 

manner  that  the  radius  is  less  than  it/2.  The  interior  of 

a circle  of  radius  r is  the  set  of  all  points  at  a distance 

2 

less  than  r from  the  center.  A disk  on  S is  a circle  of 
radius<  tt / 2 and  its  interior;  the  center  of  the  disk  is 
the  center  of  the  circle  defining  it.  (Circles  and  their 
interiors  are  not  convex  sets  if  their  radii  are  greater 
than  tt  / 2 . ) 

A . 2 Some  Topology  on 

2 

Given  any  set  A on  S a point  x is  called  an 
interior  point  of  A if  x is  the  center  of  some  disk  which 
is  entirely  contained  by  A.  A point  x is  called  an 
exterior  point  of  A if  x is  the  center  of  some  disk  which 
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is  entirely  devoid  of  points  of  A.  A point  which  is  neither 

an  interior  or  exterior  point  of  A is  called  a boundary 

point . We  denote  by  int  A , ext  A , and  bd  A , respectively, 

the  set  of  all  interior  points  of  A,  the  set  of  all  exterior 

points  of  A,  and  the  set  of  all  boundary  points  of  A, 

called  the  interior , exterior , and  boundary  of  A.  A set  K 

is  closed  if  bd  K C K (K.  contains  its  boundary  points). 

A set  is  open  if  int  K ~D  K.  (every  point  of  K is  an  interior 

2 

point  of  K) . It  easily  follows  that  for  any  set  R in  S , 
int  K = K\bd  KH  the  set  of  all  points  of  K which  do  not 
belong  to  bd  K. 

A set  K in  any  metric  space  is  called  compact 

if  for  each  sequence  of  points  {x^}  in  K some  subsequence 

{x^  } converges  in  the  metric  to  a point  x in  K (that  is, 
k 2 

lim  p(x,x  ) = 0) . Since  S is  itself  a compact  subset  of 
k+oo  nk  2 

a metric  space  it  follows  that  any  closed  subset  K of  is 

compact.  To  see  this,  assume  {x  } C K and  that  lim  p(x,x  ) 

2 n k^”  nk 

= 0 for  some  xeSz  and  some  subsequence  {x  } of  {x  }. 

nk  n 

Then  if  x/K  certainly  x is  a boundary  point  of  K:  every  disk 

with  x as  center  contains  points  in  K (namely,  the  points 

x for  large  enough  k svich  that  p(x,x  ) < radius  of  the 
uk  k 

disk),  and  a point  not  in  K.  (namely,  x itself).  Since  K 

is  closed,  bd  KC  K or  x e K , a contradiction.  Hence,  x c K 

and  {x  } converges  to  a point  in  K proving  that  K is 
nk 

compact . 


A well-known  result  of  topology  is  that  if  f is 


- 
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any  real-valued  continuous  function  defined  on  a compact 
subset  K of  a metric  space,  there  exist  points  and  in 

K such  that  f(k-^)  is  the  absolute  minimum  of  f on  K and 
f(k2>  is  the  absolute  maximum. 

Definition  A. 2 A (topological)  separation  of  a subset  A 

of  a metric  space  is  a pair  of  disjoint  open  sets  U,  V 

such  that  U and  V both  intersect  A and  AC  U U V.  A set 

A is  connected  if  it  has  no  such  separation.  (Thus,  A is 

connected  if  and  only  if  for  each  pair  of  disjoint  open 

sets  U,  V such  that  AC  U U V,  then  either  AC  U or  A C V.) 

2 

Common  examples  of  connected  sets  in  R are  points , 
segments,  arcs,  continuous  curves,  graphs  of  continuous 
functions,  and  convex  sets.  To  see  that  a segment  is  con- 
nected, let  its  endpoints  be  p and  q,  and  let  the  segment 
be  parametrized  by  { A | 0£A_<  length  of  pqHB},  with  p^tpq 
such  that  Pq  = p and  p,,  = q.  If  (U,V)  is  a separation  of 

"pq  then  without  loss  of  generality,  one  can  assume  p.  e U, 

A1 

p^  e V,  with  A ^ < A2.  Let  Aq  be  the  least  upper  bound  of 

the  set  A < A0  such  that  p e U.  Then  p.  i V,  or  else 
- L K A 0 
there  exists  a disk  D centered  at  p contained  in  V;  by 

0 

definition  of  An  a sequence  {A  } -+  A„  exists  with  p e U, 
u n u a 

n 

and  hence  lim  pA  = pA  so  pA  e DC  V for  all  sufficiently 
n-*-°°  n 0 n 

large  n,  contradicting  the  disjointness  of  U and  V.  Hence 
p^  e U.  Since  Aq  <_  A2  and  pA  e V one  must  have  Aq  f 
and  pA  e V for  ^ Let  t>e  a sequence  con- 

verging to  A _ such  that  Aq  < A^  < \„.  Then  lim  pA,  = pA 

n->°°  n 0 
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and  in  the  same  manner  that  e V was  proved  one  obtains 

n 

the  contradiction  p^ , e U.  Therefore,  no  separation  of  pq 

n 

exists  and  pq  is  connected.  (The  definition  of  convexity 
can  now  be  shown  to  imply  the  connectedness  of  any  convex 
set . ) 

The  type  argument  used  above  to  prove  the  connected- 

2 

ness  of  a line  segment  in  R will  occur  later  in  a slightly 
different  context.  First  an  obvious  property  of  connected 
sets  is  established. 

2 

Lemma  A , 1 If  A is  any  connected  set  in  R which 

contains  two  points  on  opposite  sides  of  a line  L,  then  L 

passes  through  a point  of  A. 

Proof : The  two  sides  of  L form  two  open  sets  U and 

V such  that  U U V U L = R2  . Hence  ACUUVUL.  If 

A H L = 4>  then  AC  U U V,  with  A meeting  both  U and  V. 

Hence  (U,V)  separates  A,  denying  the  connectedness  of  A. 

Recall  that  if  f is  any  real  valued  function  of  a 

real  variable,  the  graph  of  f,  denoted  graph  f,  is  the 
2 

subset  of  R consisting  of  all  points  with  coordinates 
(£»f(57)i  where  £ is  real  and  varies  on  the  domain  of  f. 

It  will  be  convenient  to  augment  this  set  in  the  case  of  a 
discontinuous  function. 

Definition  A . 3 Let  f be  any  bounded  real-valued  function 
of  a real  variable  defined  on  an  interval  a £ £ <_  8. 

The  augmented  graph  of  f is  the  set  graph  f y T,  denoted 
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by  grapn*f,  where  T consists  of  all  line  segments  joining 
(y  , f (Y ) ) with  (Y  ,Y1)  and  (V  ,Y  • where  a < y £ S and 
Y-i  = iim  sup,  f(f)  , Yo  = Iim  inf,  f(£).  (Recall  that 

c,  Y ^ '-3i 

lim  sup^^^f(C)  and  iim  inf,^f(0  are  respectively  the 
maximal  and  minimal  values  to  which  sequences  f(C  ) con- 
verge for  converging  to  y • ) Let  T denote  the 

segment  (perhaps  degenerate)  joining  the  three  collinear 
points  (y  f f (Y ) ) . (Y>Y^)  and  (y,Y2)  f°r  each  a <_  y £ 3. 


Lemma  A . 2 The  augmented  graph  in  R"  (graph*f) 
of  any  bounded  real-valued  function  f of  a real  variable 
is  connected. 

Proof : Using  the  same  notation  as  above,  for  each 

(y,6)  £ graph*f  then  a < y <_  3 and  (y,<$)  lies  on  T^  . (If 

f is  continuous  at  y,  T = the  single  point  (y,f(Y))-) 

Suppose  graph*f  is  not  connected.  Then  graph* f has  a 

2 

separation  (U,V)  with  U,V  disjoint  open  subsets  of  R , 
graph*f C U VJ  V,  and  graph*f  meets  both  U and  V.  Since 

T is  connected  for  ail  y then  either  T C U or  T,  C V 

Y Y Y 

for  ail  a <_  y _<  8 Since  all  points  of  graph*f  lie  on 

some  T then  for  certain  values  y,  and  Yo.  T CZ  U and 
Y 1 2 Y i 


T C.  V;  one  can  assume  that  < Y2-  Let  Yq  be  the  least 

upper  bound  of  all  y < such  that  T C.  U.  Then  T C U, 

- 2 Y Y0 

for  otherwise,  T C.  V,  and  bv  definition  of  y~  there  exists 
yU  u 

a sequence  {Yn}  converging  to  y such  that  Yn  £ Yq  and 

T C U.  Let  (y  ,6  ) k T ; since  {6  } is  bounded  there 
Y n n Y n 

n n 

exists  a subsequence  of  { tS n } (writing  k in  place  of  n^) 
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such  that  {6^}  converges  to  some  6^.  Hence  lim  inff  ^ f(£) 


0 

6n  < iim  supr^Y  f(0  and  (r_,6n)  e T . But  if  T d V 


- °0  - 


'0 


0’u0/  Y0 


^0 


there  is  a disk  D in  entered  at  (yq»6q)  and  contained  in 
V.  Since  lim  (yk,6k)  = (yq,60)  one  has  (yk,6kJ  e D C.  V tor 

k-+oo 

all  sufficiently  large  k.  But  (Yk>$k)  e T C U,  contra- 


Yk 


dieting  the  disiointness  of  U and  V . Therefore,  T C U 

T0 

and  since  C V,  Yq  ^ 12’  or  Y1  < Yq  < Y2’  anc*  f°r 
Yq  < Y £ Y 2 . T C.  V by  definition  of  Yq-  Now  let  (y^} 
be  a sequence  converging  to  Yq  such  that  y^  > Yq  for  all  n. 
Thus,  , C V.  Let  (Y^,5^)  f ; in  the  same  manner  as 
before  it  can  be  shown  that  for  a certain  subsequence 
( (Yk  , ) of  { ( , 6^)  ) the  point  (Yk,<$k.)  belongs  to  both 

U and  V for  all  k sufficiently  large,  in  contradiction. 

This  final  contradiction  proves  that  there  is  no  separation 
of  graph*f,  and  graph*f  is  connected. 

Observe  that  if  f is  continuous,  graph*f  = graph  f, 
and  one  obtains  the  well-known  result: 

2 

Corollary  A. 2 . 1 The  graph  in  R"  of  a continuous 
real  function  is  connected. 


A . 3 Support  Lines 

2 

Suppose  C is  a convex  set  on  S and  L is  a line 
such  that  C lies  entirely  on  L or  on  one  side  of  L (that 
is,  CC  L U H^,  where  is  one  of  the  half-spaces  deter- 
mined by  L)  and  L intersects  C at  least  at  one  point. 

Then  L is  called  a line  of  support  for  C.  Two  fundamental 
lemmas  concerning  support  lines  are  now  stated  and  proved. 


I 
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2 

Lemma  A . 3 Given  any  convex  closed  set  CC  S 
and  a point  p in  its  boundary,  there  exists  a line  of  sup- 
port for  C passing  through  p. 

Proof : If  C is  zero-  or  one-dimensional  the  claim 

is  obvious  since  CC  L for  some  line  L.  Hence  assume  C 
is  two-dimensional  and  p c bd  C.  Choose  any  line  Lq  through 
p.  If  Lq  is  not  already  a line  of  support  there  exist  two 
points  of  C on  opposite  sides  of  Ly,  q and  s,  not  on  Lq , 
and  thus  a point  r on  Lq  such  that  (qrs)  ; by  convexity  of  C, 
r e C.  Thus,  the  set  of  rays  px  for  which  x e CH  H where 
H is  the  opposite  side  of  Lq  as  s,  is  nonempty.  Accordingly, 
set 

8q  = sup{0  = •£.  xpr)  , 0 < Bq  £ ”, 

for  x e CO  H.  It  is  claimed  that  if  ray  px^  is  determined 
in  Lq  U H such  that -^Xypr  = 9q  then  line  px^  = L is  a line 
of  support  of  C through  p. 


Figure  A . 1 . 


If  L is  not  a line  of  support  there  exist 
opposite  sides  of  L,  and  thus  a point  y e 
side  of  L where  ray  pr ' opposes  pr  on  Lq . 


points  of  C on 
C lies  on  the  r'- 
If  y e H then 
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(pr'  py  pXg ) follows  from  properties  of  betweenness  (Kay 
1969)  and  0 = .$  ypr  > 4xgPr  = 0 g , contradicting  8q  1 


L 

Figure  A. 2 . 

(definition  of  Gq)  . Hence  y e Lg  U H'  x^here  H'  is  the  Xg  ' - 
side  of  Lq . But  y e Lg  implies  p is  interior  to  Ayqs  and 


Figure  A. 3. 

since  C is  convex  with  y,  q,  and  s in  C.Ayqs  and  its  interior 
belong  to  C.  Hence,  p c int  C,  a contradiction  (an  interior 
point  cannot  be  a boundary  point).  Hence  y e:  H ' and  the 
ray  py ' opposing  py  lies  on  the  H side  of  Lg . Since  (ypx') 
holds  for  some  point  x'  e C then  x'  e C O H and  Gg  > 6 ' 

=-£  x'pr.  By  definition  of  0g  there  exists  x t C f\  H such 
that  0'  < 0 < 0g  where  0 --4*xpr.  Then  it  follows  that 
(Kay  1969,  p.  63)  (px  pr'  py)  and  that  p t int  C,  a 


contradiction.  Therefore,  L is  a line  of  support  through 


Lemma  A. 4 Let  K be  a closed  subset  of  S , with 

x t K.  There  exists  p e bd  K such  that  the  perpendicular 

L to  line  at  p is  a line  of  support  of  K. 

Proof Define  the  function  f on  K by  setting 
f(k)  = p(x,k),  for  each  k c K.  It  is  well  known  that  the 
metric  itself  is  continuous  in  its  own  topology  (that 
defined  by  disks),  so  it  follows  that  f is  a continuous, 
real-valued  function  defined  on  a compact  set.  By  an 
earlier  comment  (compactness,  Section  A. 2)  f takes  on  its 
absolute  minimum  at  some  point  p e K.  Since  p(p,x)  p(k,x) 
for  all  k e K it  follows  that  p is  not  interior  to  K. 

Hence  p e bd  K.  Consider  the  perpendicular  L to  xp  at  p. 

If  L were  not  a line  of  support  there  would  exist  a point 
y e K on  the  x-side  of  L (y  \f.  L)  , and  hence  (Figure  A. 5) 

-£xpy  < tt  / 2 . Now,  x / py , for  otherwise  y e xp  which 

means  (since  x ^ y)  either  (yxp)  or  x e py  C K by  convexity 
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Figure  A . 5 . 


of  K,  a contradiction , or  (xyp)  which  would  imply  p(x,y)  < 
p(x,p),  a contradiction.  Take  first  the  case  p(x,p)  < tt/2. 
Locate  z>  the  foot  of  xj  on  line  py.  Then  since  p(x,z)  <_ 
p(x,p)  < tt/2  one  has  p(x,z)  < tt/2  and  it  follows  (Kay 
1969,  p.  114)  that  p(x,z)  < p(x,p).  Hence  z / K or  else 
the  definition  of  p e K is  contradicted.  Thus  y f z and 
p(x,z)  < p(x,y).  Let  m be  the  midpoint  of  segment  py; 
then  p(p,m)  = p(m,y)  < tt/2  since  p(p,y)  = it  implies  y r.  px. 

By  Kay  (1969,  p.  52)  either  (zpm) , z = p,  (mzp) , m = z, 

(mzy) , z = y,  or  (myz)  holds.  The  cases  z = p and  z = y are 

ruled  out  since  z / K.  Also,  since  p,  y c K and  K is  convex 

one  cannot  have  z e py;  therefore  (mzp),  m = z , and  (mzy) 
cannot  hold.  This  leaves  (zpm)  or  (myz).  If  (zpm)  then 
since  the  sides  of  triangle  xpz  are  each  < tt/2  the  Exterior 

Angle  Theorem  (Kay  1969)  gives*£mpx  > xzp  = tt/2.  But 

<mpx  =£ypx  < tt/2  (.since  y and  x are  on  the  same  side  of  L)  , 
which  is  a contradiction.  This  leaves  the  case  (myz).  Here, 
p(x,z)  < n/2  and  p(x,p)  < 7t/2  imply  p(x,m)  < tt/2  and,  in 
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turn,  p(x,y)  < rr / 2 by  Kay  (1969,  p.  98).  Therefore,  the 
Exterior  Angle  Theorem  produces  £xyp  =£xym  > xzy  - t\/2 
> £ xpy  and  by  Theorem  31.5  in  Kay  (1969),  p(x,y)  < p(x,p), 
a contradiction  since  y e K.  This  completes  the  case 
p(x,p)  < tt/2.  If  p(x,p)  > tt / 2 then  the  situation  is  as 
pictured  in  the  figure  below. 


S 

Figure  A. 6. 

Here,  choose  x'  on  with  p(x',p)  = tt / 2 and  z on  L such 
that  (x'yz).  Then  p(x,y)  < p(x,z)  + p(x',x)  £ p(x,p),  a 
contradiction.  This  is  seen  in  the  following: 
p(x,y)  _<  p(x,x')  + p(x'  ,y) 

< p ( x , x ' ) + p ( x ' , z ) 

= p(x.x')  + > (x ' ,p) 

= P(x,p) 

Therefore,  no  such  point  y c K exists  and  L is  a line  of 
support . 

2 

Lemma  A. 5 Given  a closed  convex  set  K of  S con- 
sisting of  more  than  one  point  and  a point  x / K such  that 
no  line  through  x contains  K,  there  exist  exactly  two 


1 


distinct  rays  xb-^  and  xb2  with  b^  e bd  K,  k = 1,2,  such 
that  the  lines  *xb^  are  support  lines  , and  the  ray  5cu  for 
each  u c K either  coincides  with  xb^ , or  xb^ , or  lies 
between  xb^  and  xb^ . 

Proof : Since  no  line  through  K.  contains  x and  K 
is  not  a point,  there  is  a segment  y ~y~>  Cl  K such  that 
z e yiy2 , and  y^  and  y^  lie  on  opposite  sides  of  line  xz . 
Define  the  bounded  nonempty  sets: 

S^.  ={(1,  ^tr.xu^j  u,  lK,  u^  lies  on  the  y^-side  of  xz , 
k = 1,2} 

Let  O*  < n be  the  least  upper  bound  of  , k = 1,2.  Choos- 

ing  a sequence  {b^}  such  that  bn  c xu^i'  K for  eacb  k,  a 

subsequence  (b^'l  converges  to  some  point  b,  e K (since  K 
m k 

is  compact)  such  that  9,  =^Czxb,  , k = 1,2. 


Figure  A . 7 . 


Observe  that  b^  i bd  K,  for,  if  b^  { bd  K then  b^ 
is  an  interior  point  of  K,  in  which  case  a disk  D about  b. 
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lies  in  K;  but  then  there  exist  u c D H V such  that  ray 
xb\.  lies  between  xu  and  xz  and  u is  on  the  y^-side  of  xz 
with  0 = * zxu  ^^zxb^  d-^b^xu  > -4  zxb^  = 0^,  contradicting 
the  fact  that  6^  is  an  upper  bound  for  . 

Now,  if  xb^.  is  not  a support  .line  of  K there  exists 
points  u^  , u2  of  K on  opposite  sides  of  line  xbj.  , and  the 
segment  uJiT,  intersects  xb^.  at  a point  u^  ; if  u^  / ray 
xb^  or  u^  = x then  x e conv  { u^  ,u2,z}  C K contradicting 
x / K.  Hence  u^  c xb^.  and  u^  ^ x,  so  that  a point  u e 
~2  C\  K such  that  ray  xb^  lies  between  xtT  and  "xz  and  u 

4— > * 

is  on  the  y^  side  of  xz  with  0 = -41  zxu  >-izxbj<.  = 6^ 

(as  before),  in  contradiction.  Therefore,  xb,  is  a support 


line  of  K 


Finally,  if  u is  any  point  of  K then  u lies  both 


on  the  b2  side  of  xbp  and  the  b^  side  of  xb2-  That  is, 


u e [H(b^x,b2)  O H(b2x,b^)]  U xb'  U xb  2 and  u lies  on  xb^ 

xb*,  , or  is  an  interior  point  of  the  angle  b^xb2>  Hence 

either  xu  coincides  with  xb^  or  xb*  , or  jcu  lies  between 
--  ^ 

xb  ^ and  xb  9 . 

(Note  that,  in  particular,  z e K.  so  (xb^  xz  xb2) 
or<9.zxb^  +^.zxb2  =4b^xb2  £ 7t  and  6,  + 0O  £ n.) 


Lemma  A . 6 Given  a closed  convex  set  V of  S , with 
points  x / V and  b t bd  V such  that  tc?  is  a support  line 
to  V at  b.  There  exists  a point  x'  on  bd  V and  u f xb  such 
that  ux'  is  a support  line  of  V at  x'  and  p(u,x')  = p(x,u). 
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Proof : Let  g be  defined  as  follows : For  each 

u t "xE  let  u'  denote  that  point  on  L D bd  V nearest  u, 

r u 

where  Lu  denotes  the  unique  line  of  support  of  V from  u 

different  from  the  support  line  xb  (Lemma  A. 5).  Then  take 

g(A)=  g(u.)  = p(uA,up  - p(x,uv),  0 <_  A £ 0,  where  uA 

denotes  a parametrization  of  segment  xb  such  that  p(b,uv)  = A, 

with  un  = b and  uQ  = x.  Note  that 

u P 

g (0)  = p (b  , b ' ) - p ( x , b ) = -p(x,b)  < 0,  g (8 ) = p(x,x')  - p(x,x) 
= p (x,x ' ) > 0 . 

2 

The  augmented  graph  of  g in  R (coordinated  by 
(£,n))  is  connected  by  Lemma  A. 2,  with  the  points  (0,g(0)) 
and  (8,g(3))  on  opposite  sides  of  the  line  n = 0 (the  -axis); 
hence  (y,0)  E-  graph*g  for  some  y.0£y£8.  ifgis 
continuous  at  y then  iim  inf,^g(A)  = g(Y)  = lim  supA 
g ( A ) and  (y,0)  = (Y,g(y)).  In  this  case  one  has  g(Y>  = 0, 
or,  p ( u^  ,u')  = p(x,u  ) . 

However,  g may  be  discontinuous  at  Y:  so  it  is 
desired  to  show  that  g has  a simple  jump  discontinuity  at 
any  such  point.  From  the  nature  of  the  boundary  of  a convex 
set  it  is  clear  that  for  points  uA  with  A < y the  corres- 
ponding point  uj  on  bd  V precedes  or  coincides  with  u*  , 
taking  a counterclockwise  orientation  (parametrization) 
of  bd  V.  For,  if  u’  follows  u,1  on  bd  V (Figure  A. 8)  then 

bu!  meets  line  u u’  at  a point  v such  that  either  v = u! 

A V Y A 

* — v 

or  u!  is  between  v and  b,  since  u u*  is  a line  of  support 
A y Y 

of  V and  the  points  b and  u[  cannot  lie  on  opposite  sides 
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I 

| 

| 

| 

| 

I 


of  n^u^.  By  the  Crossbar  Principle  (Kay  1969)  and  related 
ideas,  the  segment  meets  u,|b  at  some  point  w interior 

to  both  segments,  and  hence  V contains  points  u^'  and  b on 
opposite  sides  of  line  u^u^,  contradicting  the  fact  that 


Figure  A. 8. 


Now  let  y be  a point  of  discontinuity  of  g. 

Let  A>y~  (that  is,  A approaches  y from  values  less  than  y) . 

Since  bd  V is  compact  some  sequence  {A^}  of  A's  is  such 

that  {u!  ) converges  to  a point  v,  bd  V.  If  lim  u ' / v, 
n 1 A-*°° 

then  another  sequence  of  A's  is  such  that  fu^  } con- 


n 


verges  to  v9  / Vj  on  bd  V.  But  the  line  u^  uj  converges 

n n . 

to  the  line  u v-^ , and  this  Line  must  be  a line  of  support 
to  V from  u . Since  there  are  exactly  two  such  support  lines 
(Lemma  A.  5),  u v-,  must  be  either  L or  xb  ; it  is  obvious 

V X U 


that  u v,  = L 
Y 1 u 


^ — ■*» 


Similarly,  u v0  = L and  hence  both  v, 

7 ’ Y 2 u 1 

t * — ^ y 

and  v0  lie  on  Lu  = u^u'  with  and  v coinciding  with  or 
preceding  u^  on  bd  V.  But  by  definition  of  u^ , 


I I 
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P(uy  ,vk)  > P(uY  ,uj)  , k = 1.2 
and  if  inequality  prevails  then  vk  would  follow  u^  on  bd  V, 
which  is  impossible. 


Hence 


P(UY.V  = p(uy 


with  v,  e L . Hence  v,  — u'  for  k — 1,2.  This  proves 
k u k y 

Y 

that  lim  u*=  u'  and  hence 

^ ^ A y 

lim_  g(A)  = lim_  tp(u^,u^)  - p(x,u^)]  = p(u^,u^)  - p(x,u^)  = g(Y) 
A+y  A-*-y 

By  an  argument  exactly  similar  it  can  be  shown  that 

if  u"  is  that  point  on  L A bd  V farthest  from  u then 
Y uy  Y 

lim  u*  = uy' 

A-Y 


lim,  g ( A ) = p(u  , u";  - p(x,u  )=r  a 


, + ' Y Y 

A->y 

Since  p(u  ,uY)  > pCu^.u^l  we  have 

p(uY , u”)  - p(x,uA)>  g ( Y ) 


lim , g ( A ) = ot  _>  g ( Y ) • 

A-"Y  Y 

Hence  if  g is  discontinuous  at  y then 

lim  g ( A ) = a > g(Y)  = lim_  g(A)  . 

A-*-y  ' a-»y 

Thus  it  follows  that 

lim  infA+yg(A)  = mi n { g( y ) , <\  1 = 

lim  sup  g(A)  = max{g(Y),Uy}  = <*Y 

and  (y.O)  c graph*g  implies  (y,0)  lies  on  the  segment 

determined  by  the  points  (y , g( Y) ) , (Y.oO  and  g(Y)  £ 0 1 ay • 


That  is, 


p(uY,u^)  - p(x,uy)  < 0 < p(u^ ,u")  - p(x,uy) 
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p(uY,u^)  <_  p(x,uY)  £ p(u  ,u^'). 

Since  pCu^.v)  varies  continuously  from  p(u^,u^)  to  p(u^,u”) 
as  v varies  on  the  segment  u^u" , the  ordinary  Intermediate 
Value  Theorem  guarantees  a point  v on  u^u"d  bd  V such  that 

p(uy,v)  = p(x.U^) . 

Thus  in  all  cases  it  has  been  proven  that  there  exists  a 
point  x'  on  bd  V and  u such  that  ux'  is  a support  line  of 
V at  x'  and  p(u.x')  = p(x,u). 


Definition  A. 4 The  Caratheodory  number  c is  the  least 

integer  for  which  it  is  true  that  for  any  set  A and  any 

point  x e conv  A there  exist  c points  a-^ a£  of  A 

such  that  x e convla, , ....  a 1 . 

1 c 

Lemma  A . 7 The  number  c (Caratheodory)  is  equal  to 
3 in  spherical  space,  as  it  is  on  the  plane. 


Proof : It  is  clear  that  c > 3 by  considering  the 

case  A = {a^,  a^l , the  vertices  of  a nondegenerate 

spherical  triangle.  Hence  it  suffices  to  show  that  c 3 . 
2 

Let  A d S and  suppose  x r.  conv  A;  also,  suppose  x / A, 
for  otherwise  one  may  simply  take  a^  = a9  = a^  = x and  the 
result  is  trivial.  Assume  first  that  A is  contained  by 
some  closed  hemisphere.  Consider  the  two  cases:  x c bd 

conv  A and  x c int  conv  A.  If  x c bd  conv  A let  L be  a 
line  of  support  of  conv  A at  x and  H the  A-side  of  L if 
A^.L;  here  it  follows  that  x t conviap  a9f  for  a pair  of 
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points  a^,  c A (and  we  take  a^  = a2  so  that  x e conv 
[ap  a2 , a^})-  For,  if  no  point  of  L lies  in  A then  H is 
a convex  set  containing  A but  not  x,  which  implies  the  con- 
tradiction x e conv  AC  conv  H = H,  since  x 4 H.  Thus, 
there  must  exist  a e A H L.  Consider  the  two  opposing  rays 
on  L with  origin  x,  say  and  R2 , with  a >.  R^.  Now  A C\ 

R2  ^ <P  for  otherwise  H KJ  RL\  ( x}  is  a convex  set  containing 
A but  not  x.  Let  a^  t A O R^,  k = 1,  2,  be  such  that  f (x.aj  ) 
is  minimal  (assuming  for  the  moment  that  A is  closed;  the 
argument  for  the  case  when  A is  not  closed  is  very  similar) . 

If  p(x,a-^)  + p(x,a0)  71  then  x 4 convfa-^,  a2)  and  convla^  , a0 

contains  all  the  points  of  A C\  L (by  the  minimal  property 
of  p(x,a^),  k = 1,2);  then  one  would  have  the  convex  set 
HU  convla-^,  a2)  containing  A but  not  x,  a contradiction. 

Hence  p(x,a^)  + p(x,a7)  < it  and  it  follows  that  x c a ^a7  = 
conv(a^,  a2^ ■ This  completes  the  case  x e bd  conv  A. 

If  x e int  conv  A then  let  D be  a disk  centered 
at  x with  DC  conv  A.  There  must  exist  a^  c A such  that 
a^  lies  in  the  open  hemisphere  H whose  closure  contains  A, 
for  otherwise  A would  be  contained  in  the  line  L which 
determines  H and  conv  ACL  or  D <£.  conv  A,  a contradiction. 

Hence  a^  c AO  H and  we  consider  line  xa^  and  the  two 

^ ^ 

sides  H^  and  H2  of  xa^  . Certainly  A f\  H^  / for  k = 1,2, 
or  else  conv  AC  conv(A  H H^+j)  and  D conv  A.  The  clo- 
sure of  A f\  H,  (the  set  A f\  H^  and  all  its  boundary  points, 
denoted  cl(AO  H^))  is  compact,  so  there  can  be  found 
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a£  e cl(A  A H^)  and  a^  e cl(AO  H^)  such  that  4 a^xa2  an<^ 
^a^xa^  are  maximal  (Figure  A.  9).  If  4a^xa2  +-4-a^xa3  1 77 
then  all  of  A would  lie  on  one  side  of  a line  through  x 
and  we  would  have  D <£.  conv  A.  Hence  4a^xa2  + 4aixa3  > 71 
and  there  accordingly  exist  a.2  e A C\  H^,  a^  c A C\  such 

that  4 3-^2  + 4-a^xa^  > tt  (Figure  A. 9).  It  follows  since 


Figure  A. 9. 

a^  £ H and  a^,  a^  e HO  L,  that  x is  an  interior  point  of 
triangle  a-^a^a^.  This  implies  that  x e convia^,  a^,  33)  , 
finishing  the  case  when  A lies  in  a closed  hemisphere. 

Finally,  if  A does  not  belong  to  a closed  hemisphere  then 
let  H be  a closed  hemisphere  which  contains  x,  and  define 
the  set  A'  = A f\  H.  Then 

conv  A'  = conv(A  f\  H)  = conv  A O conv  H = (conv  A) A H 
and  since  x belongs  to  both  conv  A and  H we  have  x c conv  A' 
where  A'  is  contained  by  a closed  hemisphere.  By  the  pre- 
ceding case,  x e convta^,  a^,  a where  a^  , a0,  a^  belong 

to  A',  and  therefore,  also  to  A. 
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A . 4 Convex  Polytopes  on  S 

A convex  polytope  is  the  convex  hull  of  a finite 

set  of  points.  As  in  classical  convexity,  every  convex 
2 

polytope  in  S is  compact.  The  question  of  when  conv  A is 

2 

compact  is  nontrivial  even  though  on  S any  closed  set  is 
compact.  Just  as  in  the  plane,  there  exist  closed  sets  A 
such  that  conv  A is  not  closed:  Take  a pair  of  antipodal 

points  a,  a'  and  a circle  C through  a'  of  radius  less  than 
tt / 2 (Figure  A.  10).  Then  take  A = aU  C;  the  convex  hull 
of  A is  clearly  the  open  hemisphere  containing  C,  together 
with  a and  a'--not  a closed  set. 

a ' 

2 

Figure  A . 10 . 

It  is  now  shown  that  if  A is  finite  then  conv  A 
is  closed  (therefore  compact) . 

Lemma  A. 8 Any  convex  polytope  on  is  compact. 

Proof : Let  A = {a^,  ..,  anl , and  consider  a boundary 

point  x of  conv  A;  it  is  shown  that  x e conv  A,  proving 
that  conv  A contains  its  boundary  points.  For  each  integer 


1 
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k = 1,2,  ...  let  be  the  disk  centered  at  x having  radius 

1/k.  Then  D,  contains  a point  x,  e conv  A,  and  lim  x,  = x 

, , k^ 

(since  lim  pCx.x^)  = 0).  Since  S has  Caratheodory  number 
k-*-°° 

3 (Lemma  A. 7)  there  exist  for  each  k a set  of  3 points 
a,  , a0  , and  a0  from  a, , ....  a such  that  x,  e conv 

Y V \ 1 k 

{a1  , a9  , a9  } . Since  A is  finite  there  is  a subsequence 

lk  k \ 

{ k ' } of  k's  for  which  a-,  is  constantly  equal  to  a in  A. 

V u 

Hence  x.  , c convla  , a,,  , aQ  }.  Again,  a subsequence  {k"} 

k u 2k’  3k' 

of  { k ' } exists  such  that  a0  = a e A and,  finally,  a 

Zk"  v 

subsequence  { i } of  {k"}  such  that  aQ  = a . Hence,  a 

j ^ w 

subsequence  1 i } of  { k > exists  such  that  x.  c convla  , a , a }, 

1 U V w 

with  lim  x.  = x.  Since  conv  la  . a , a } is  either  a 

. 1 U V W 

j_-+-oo 

point,  a segment,  a line,  or  a nondegenerate  triangle 

and  interior  ("therefore  closed),  it  follows  that  x e 

convla  , a , a } CL  conv  A,  as  was  to  have  been  proved, 
u ’ v w r 


APPENDIX  B 

DATA  FOR  EXAMPLE  PROBLEMS 


n , Latitude 

Number  point  Longitude  Weight 

in  Degrees 
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n , Latitude 

Number  eman  Longitude  Weight 
romt  in  Degrees 
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Latitude 

Number  Demand  Point  Longitude  Weight 

in  Degrees 


D16 

See 

Figure 

5.2 


(-12,48) 

(65,75) 

(15,-20) 

(25,-115) 

(-30,175) 

(-70,-110) 


P, -Norway 
Station 

(-70,0) 

1 

P„-Amundsen- 
Scott  Stn. 

P~-Palmer 

Station 

(-90,0) 

(-65,-65) 

1 

1 

Lj 

/ 

0°  Meridian 

P, -Byrd 
Station 

D17  Pc-Little 

America 

(-80,-120) 

(-78,-160) 

1 ^ 

vj 

} 

^V7 

1 

Y 4# 

P, -Port 
Martin 

(-66,140) 

1 

\ — 

P5 

\js/ 

P ^ -Maws on 

(-68,62) 

1 

L80° 

Pg-Mirnyy 

(-66,94) 

1 

Pq-American 

Highland 

(-70,75) 

1 
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Number  D18: 


City 

Latitude 

Longitude 

Weight 

City 

Latitude 

Longitude 

Weight 

Minsk 

(54,28) 

0.012 

Zhdanov 

(47,38) 

0.007 

Lvov 

(50,24) 

0.010 

Stalino 

(48,37) 

0.016 

Kishinev 

(47,28) 

0.005 

Makeyevka 

(48,38) 

0.008 

Odessa 

(46,31) 

0.014 

Gorlovka 

(48,38) 

0.007 

Nikolaev 

(47,32) 

0.005 

Taganrog 

(47,39) 

0.005 

Kherson 

(4/, 33) 

0.004 

Krasnodar 

(45  39) 

0.007 

Sevastapol 

(45,34) 

0.004 

Rostov 

(47,40) 

0.014 

Simferopol 

(45,34) 

0.U04 

Shakhty 

(48,40) 

0.004 

Krivol  Rog 

(46,34) 

0.009 

Kadiyevka 

(49,39) 

0.004 

Dneprodzerhinsk 

(48,35) 

0.004 

Kharkhov 

(41,37) 

0.021 

Dnepropetrovsk 

(48,35) 

0.015 

Kiev 

(50,31) 

0.026 

Zaporozhe 

(48,36) 

0.010 

Gomel 

(52,31) 

0.004 

Number  D19 : 

Demand 

Point 

City 

Latitude 

Longitude 

Weight 

1 

Montgomery,  Alabama 

(32°23 ' 

, -86°17 ' ) 

62,650 

2 

Juneau,  Alaska 

(58°25 ' 

, - 134°30 ' ) 

1,163 

3 

Phoenix,  Arizona 

(33°30 ' 

, -112°00 ' ) 

22,288 

4 

Little  Rock,  Arkansas 

(34°42 ' 

, -92°16 ’ ) 

20,916 

5 

Sacramento.  California 

(38°35 ' 

, -121°30 ' ) 

394,139 

6 

Denver,  Colorado 

(39°44 ' 

, - 104°59 ' ) 

47,453 

7 

Hartford,  Con  "cticut 

(41°45 ' 

, - 72°40 ' ) 

74,813 

8 

Dover,  Delaware 

(39°10 ' 

, -75°30 ' ) 

15,863 

9 

Washington,  D.C. 

(38°50 ' 

1 

o 

o 

o 

238,476 

10 

Tallahassee,  Florida 

( 3U°25  ' 

, -84°1 7 ' ) 

80,791 

11 

Atlanta,  Georgia 

(33°45  ' 

, - 84°23 ' ) 

79,198 

12 

Boise,  Idaho 

(43°38 ' 

, -116°1 2 ’ ) 

9,502 

13 

Springfield,  Illinois 

( 39°46 ' 

, -89°37 ' ) 

387,961 

14 

Indianapolis,  Indiana 

( 39°45 ' 

, -86°08 ' ) 

91,294 

15 

Des  Moines,  Iowa 

(41°35 ' 

, -93°37 ' ) 

60,204 
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Demand 

City 

Latitude 

Weight 

Point 

Longitude 

16 

Topeka,  Kansas 

( 39°02 ' 

, -95°41 ' ) 

47,086 

17 

Frankfort,  Kentucky 

(38°10 ' 

, -84°55 ' ) 

41,344 

18 

Baton  Rouge,  Louisiana 

(. 30°28 ' 

, -9 1°10 ' ) 

52,978 

19 

Augusta,  Maine 

(44°19 ' 

, -69°42 ' ) 

17,631 

20 

Annapolis,  Maryland 

( 39°00 ' 

, -76°25 ’ ) 

73,834 

21 

Boston,  Massachusetts 

(42°15 ' 

, -71°07 ' ) 

153,409 

22 

Lansing,  Michigan 

(42°45 ' 

, -84°35 ' ) 

127,187 

23 

St.  Paul,  Minnesota 

(44°57 ' 

, -93°05 ' ) 

90,323 

24 

Jackson,  Mississippi 

( 32°17 ’ 

, -90°10 ' ) 

24,448 

25 

Jefferson  City,  Missouri 

(38°34 ' 

, -92°10 ' ) 

139,140 

26 

Helena,  Mon t ana 

(46°34 ' 

, - 112°01 ' ) 

17,322 

2/ 

Lincoln,  Nebraksa 

(40°49 ' 

, -96°43 ' ) 

35,658 

28 

Carson  City,  Nevada 

(39°10 ' 

, -119°45 ' ) 

9,207 

29 

Concord,  New  Hampshire 

(43°1(J  ' 

, - 71°30  ' ) 

11,631 

30 

Trenton,  New  Jersey 

(40°13 ' 

, - 74°46  ' ) 

184,397 

31 

Sante  Fe , New  Mexico 

( 35°10  ' 

,-l06°00') 

17,645 

32 

Albany,  New  York 

(42°40 ’ 

, -73°50') 

662 ,584 

33 

Raleigh,  North  Carolina 

( 35°45 ' 

, - 78°39  ' ) 

73,749 

34 

Bismark,  North  Dakota 

(46°48 ' 

, -100°46  ' ) 

11,646 

35 

Columbus , Ohio 

(40°00 ' 

, -83°00 ' ) 

219,330 

36 

Oklahoma  City,  Oklahoma 

(35°27 ' 

, -97°32 ' ) 

59,159 

37 

Salem,  Oregon 

(44°55 ' 

, -123°03 ' ) 

45,733 

38 

Harrisburg,  Pennsylvania 

(40°15 ' 

, -76°50 ' ) 

302,933 

39 

Providence,  Rhode  Island 

(41°50 ' 

, - 71°23 ' ) 

22,769 

40 

Columbia,  South  Carolina 

(34°00  1 

, -81°00 ' ) 

28,434 

41 

Pierre,  South  Dakota 

(44°22 ' 

, - 100°20  ' ) 

10,292 

42 

Nashville,  Tennessee 

( 36°10 ' 

, -86°48 ' ) 

68,770 

43 

Austin,  Texas 

( 30°15 ' 

, -9  /°42 ' ) 

233,041 

44 

Salt  Lake  City,  Utah 

(40°45 ' 

, -111°52 ' ) 

23,110 

45 

Montpelier,  Vermont 

(44°2U ' 

, -72°35 ' ) 

14,082 

; 

46 

Richmond,  Virginia 

( 37°35 ' 

, - 77°30  ' ) 

74,408 

Demand 

Point 

City 

Latitude 

Longitude 

Weight 

47 

Olympia , 

Washington 

(47°02 1 , -122°52 ' ) 

53,472 

48 

Charleston,  West  Virginia  ( 38°20 -8i°35 ' ) 

23,240 

49 

Madison , 

Wisconsin 

(43°05 ' , -89°23 ' ) 

86,544 

50 

Cheyenne 

, Wyoming 

(41° 10 ' , -104°49 ' ) 

7,489 
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